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Preface 


This two-part document contains text and figures for the papers presented 
at the Symposium on Transonic Unsteady Aerodynamics and Aeroelasticity — 1987, 
held at the NASA Langley Research Center on May 20-22, 1987. The Symposium, 
which reviewed the subject area, was the third such meeting at Langley. The 
previous workshops were held in 1980 and 1983 and were an outgrowth of a meet- 
ing held at Columbus, Ohio in 1978 to assess the state of unsteady aerodynamics 
for use in transonic flutter analysis. 

The papers were grouped in five subject areas, which may be described 
broadly as: 

(1) Transonic small disturbance (TSD) theory for complete aircraft 
configurations 

(2) Full potential and Euler equation methods 

(3) Methods for vortex and viscous flows 

(4) Aeroelastic applications 

(5) Experimental results and cascade flows 

The decade since the Columbus meeting has seen the wide acceptance of 
computational fluid dynamics methods for transonic aeroelastic analysis. In 
1978, calculations with the TSD methods for two-dimensional airfoils (espe- 
cially the NASA Ames LTRAN code and its derivatives) were well established, 
and the USAF Flight Dynamics Laboratory had initiated development of the first 
TSD code (XTRAN3) for three-dimensional wings. As demonstrated in the papers 
from the 1987 meeting contained herein, the TSD methods (the NASA Langley 
CAP-TSD code in particular) can now be applied to the aeroelastic analysis of 
complete aircraft. Methods suitable for situations in which small disturbance 
theory is inadequate are being aggressively pursued. 

Future research should follow three main paths: 

(1) Development of more exact methods using the full potential, Euler, 
and Navier-Stokes equations 

(2) Evaluation of the TSD methods by detailed comparison with more 
exact methods and experiment 

(3) Detailed pressure measurements and flutter tests on well-defined 
aeroelastic models. 


Samuel R. Bland 
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OUTLINE OF TALK 


1. Background and Objectives 

2. Unsteady Euler Equations in a Rotating Frame of Reference for Transonic Vortex 
Flows: 

- Formulation (Space Fixed and Rotating Frames of Reference) 

- Method of Solution, Local-Conical Flow, Initial and Boundary Conditions 

- Results: • Symmetric Conical Flow 

• Three-Dimensional Steady Transonic-Vortex Flow 

• Uniform Rolling in a Conical Flow 

• Rolling Oscillation in a Conical Flow 

3. Transonic Airfoil Computation Using Integral Solution of Full-Potential Eq. With 
and Without Embedded Euler Domains: 

Formulation (I.E. Solution of Full-Potential Eq., Euler Equations) 

- Method of Solution (SCSF-Scheme , IEEE-Scheme) 

- Results: SCSF and IEEE Schemes are Applied to 

NACA 0012 and NACA 64A010A Over a Wide Range of M . 

oo 

4. Concluding Remarks 
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1. Background and Objectives 


• As the normal angle of attack, normal Mach number and sweep-back angle are 
varied; complex flows develop around Delta and belta-like wings. 

These flows are characterized with the formation of large and small scale 
vortices , weak and strong shock waves , and shock induced separations . The 
experimental Miller and Wood Classification Diagram shows seven regions of this 
flow. 


• These flows become highly complex when vortex breakdown occurs in the vicinity of 
the wing or when the wing undergoes unsteady motion due to maneuvering or 
flutter. 


• The main objectives of this ongoing research research work are to develop 

efficient and reliable computational schemes which are capable of predicting the 
distributed aerodynamic characteristics of these wings in steady and unsteady 
flows over a wide range of angles of attack, sweep-back angles, Mach numbers and 
configurations. 





1 Classical 
Vortex 

2 Separation 
Bubble with 
Mo Shock 

3 No Shock/ 

No Separation 

4 Shock with no 
Separation 

5 Shock- Induced 
Separation 

6 Separation 
Bubble with 
Shock 

7 Vortex with 
Shock 


Fig. l Miller and Wood 1 Classification Diagram. 
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2. Unsteady Euler Equations in a Rotating Frame of 
Reference for Transonic-Vortex Flows 

Formulation 

• Conservation Form of Euler Equations in a Space-Fixed Frame of Reference 


*i + &!+ o 

(1) 

at ax ay az 

where 

q = [p, pu, pv, pw, pe] 1 

(2) 

2 t 

F = [pu, pu + p, puv, puw, puh] 

(3) 

F = [pv, puv, p v 2 + p, pvw, pvh] t 

(4) 

2 _,t 

<a = [pw, puw, pvw, pw + p, pwhj 

(5) 

e = p/p(y-l) + (u 2 + v 2 + w 2 )/? 

(6) 

h = e + p/p 

(7) 

• Rewriting the Equations in the Vector Form 


+ 7 . (p V) » o 

(8) 

+ v.(pVv*p!)«o 
9t 

(9) 

M + V • (p h V) - 0 

(10) 
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And Using the Substantial and Local Derivatives Relations 


Da D'a aa a' a 

Dt Dt ' at at v r ' 


(ID 


DA D'A - - aA a 'A T - T 


( 12 ) 


where V 
Frame 


- - D' a' 

V + wxr ; — - , — = Substantial and Local Derivative in the Rotating 
r nt 


• We get the Conservative Form of Euler Equations for the Relative Motion 


^ + V.(pV r )=° 


(13) 


a'( P v p ) 

— + V • [p V V + pl] = - p[wxr + 2cjxV + wx(wxr)] 

ot r r r 


(14) 


9' (p e p ) _ * _ • 

— + v • [p h r V^] = - p[V r • (wxr) + (uxr) • (wxr)] 


(15) 


where 


P V r 1 |- - ( 2 - - x 

e r = tmT) + ~- 2 l“ Kr l = e - v ‘ (wxr) 

V ? 1 _ _ 2 

h = , yP .\ + - — Lxrl = h - V • (uxr) 


r p(y-l) 2 2 


(16) 


(17) 
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Method of Solution 


• The Abstract Conservative Form of the Relative Motion in Terms of Rotating 
Coordinates are 


a'q a'E a'F a'G 

r r r r - 

+ + + = s 


at 


ax ay az 


( 18 ) 


where 


q r = Cp, pu r , pv r , pw r , pe r ] t 


(19) 


F = 
r 

F = 
r 

G = 
r 

S 


_ 9 

2 

P u + p, pU v , puw, 

PU h f 

(20) 

r 

r r r r r 

r r 


_ 9 

2 

pu v , pv + p, pvw. 

pv h T 

(21) 

r 

r r r r 

r r 


r’ 

2 

pu w , pV W , pw + p, 
r r r r r 

pw h 
p r r 

(22) 


[0, 0, p(uz + 2w w^ + «y), - p(coy + 2cov^ 


w 2 z), 


p(-v o)Z + v/ oiy + wojy 
r r 


+ oxLz )] 


2 \ -it 


(23) 


The Source Term 5 has been Written for JL = co e , w = w e 

a A 

• Eqs. ( 18) - ( 23 ) , (16) and (17) are Solved Using a Central-Difference, Finite-Volume 
Scheme Using Four-Stage Runge Kutta Time Stepping with Added Second- and Fourth- 
Order Dissipation Terms. 

• The Computer Program is a Three-Dimensional Program. 
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• Local-Conical Flow Problem 


• If the Conical Coordinates are Used to Transform the Relative Motion Equation, 
Eqs. (18)-(22), the Resulting Equations will not Represent a Conical Flow. 


® If the Conical Coordinates are Used to Transform the Absolute Motion Equations, 

Eqs. (l)-{7), the Resulting Equations will Represent a Conical Flow for the Steady 
Flow. For the Unsteady Flow, the Problem is Made "Locally Conical" if it is Solved 
at a Fixed Axial Location. 


• In the Three-Dimensional Program, Local Conical Flow Solutions are Obtained at x=l 
by Equating the Absolute Motion in the First and Third Planes: 


(pV) 1tl - (pV), - (pV r ) tll - (pV r ) = (pV r >. + u> x( Pj f t - P(tl ? 1tl ) 


(24) 


<P E) 1t l ■ ((«), * <P E r ) 1± l ' (P E r >, + < 5 r 1 + ? ( - < 25) 


where i = 2. 
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1. Initial Condition: 


For Constant Rolling Problem (w = -w e ),_the Flow Corresponds to a Uniform 
Translation Plus a Rigid Body Rotation -ijxr. 


For Rolling Oscillation (w = -w sin kt e ), the Flow Corresponds to a Uniform 
Translation Only. 


2. Boundary Conditions : 


Normal Momentum Equation is Used on the Wing Surface 
({Sr (v r • n ) = 0, = 0) 


pV . (V • vn) = — + p n • [2wxV + wxr + ojx(o)xr)^ (26) 


In the Farfield, a Uniform Translation Plus the Corresponding Rotation are 
Imposed Outside of the Bow Shock. 
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SYMMETRIC COMICAL FLOW 


This is a verification test case for the three-dimensional program which has 
been solved earlier by using a conical flow program^. Figure 2 shows the results 
for a flat plate sharp-edged delta wing at M m = '2, a = 10° and e (sweep angle) 

using a modified Joukowski transformation^ of 128x64 cells around and normal to the 
wing has been used for the whole computational region. Figure 2 (a--d) shows the 
surface pressure, cross-flow velocity, cross-flow Mach contours and static pressure 
contours. It is clear that two symmetric leading-edge vortices have been captured on 
the suction side along with a weak cross-flow shock under each vortex. The outer bow 
shock is clearly visible in the lower portions of the cross-flow Mach and static 
pressure contours. The surface pressure, Mach contours and static-pressure contours 
match those obtained by using the conical flow program for half of the computational 
region. 



a. surface pressure b. cross-flow velocity 

Fig., 2 Steady symmetric flow around a delta 
wing, M oo =2, a =10°, p=70° 


269 






ORIGINAL PAGE IS 
OF POOR QUALITY 


THREE-DIMENSIONAL TRANSONIC FLOW 


Figure 3 shows the results of a sharp-edged delta wing for M^ = 0.7, a = 10° and 
aspect ratio of 1.5 using a number of cells of 80x38x48 in the x,°° n and $ direc- 
tions; respectively. The results of Figures 3a and 3b show, from left to right, 
the surface pressure, the static pressure contours and the cross-flow velocity; each 
at the chord stations of 0.52 and 0.81. Figure 3c shows the static pressure 
contours and the cross-flow velocity at the chord station of 1.01 and Figure 3a 
shows the static pressure contours and the cross-flow velocity at the chord station 
of 1.25. At x = 0.81, comparisons of the computed surface pressure with the 

experimental data^ shows that the location of vortex core is well predicted, while 
the value of peak suction pressure under the vortex core is slightly under- 
predicted. At x = 1.25, the static pressure contours and the cross-flow velocity 
show the formation of trailing-edge vortex core and its interaction with the leading- 
edge vortex core. Our cross-flow planes are taken normal to the wing surface. 



a. X =0.52 


b. X =0.01 


Figure 3 . Three-Dimensional Transonic Flow, Standard Euler Set, Sharp-edged delta 
wing, 80X38X48 cell, M^O.7, o»15°, AIM. 5, e 2 =0.12, t^O.QOS, 

(».,b.) 1. Surface Pressure, 2. Static Pressure Contours, 3. Crossflow 

velocity 

Static pressure Contours, 2. Crossflow velocity 
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Figure 3 . Three-Dimensional Transonic Flow, Standard Euler Set, Sharp-edged delta 
wing, 80X38X48 cell, M^O.7, o*15°, AR*1.5, e 2 =0.12. e 4 =0.005. 

(a, b.) 1. Surface Pressure, 2. Static Pressure Contours, 3. Crossflow 

velocity 

(c.,d.) i. static pressure Contours, 2. Crossflow velocity 
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THREE-DIMENSIONAL LOW-SPEED FLOW 


ORIGINAL PAGE IS 
OF POOR QUALITY 


Figure 4 shows the results of a sharp-edged delta wing for M^ = 0.3, a = 20.5° 
and aspect ratio of 1 using a number of cells of 80x38x48 in the x, n> and ? 
directions; respectively. The results of Figures 4a and 4b show, from left to 
right, the surface pressure and the experimental data of Hummel^ the static-pressure- 
coefficient contours and the cross-flow velocity at the chord stations of 0.52 and 
0.81. The results of Figures 4c and,4d show, from left to right, the experimental 

static-pressure-coefficient contours^ (cross-flow planes are normal to wind 
direction), the computed static-pressure-coefficient (cross-flow planes are normal to 
wing surface) and the computed cross-flow velocity at the chord stations 1.02 and 


1.25. 


The results show that the location of leading-edge vortex core is accurately 
predicted, the suction pressure peak is accurately predicted at x = 0.52 but it is 
overpredicted at x = 0.81, the location of the trai 1 ing-edged vortex core is slightly 
off that of the experimental data. A finer grid than the one used in this example is 
expected to give a more accurate prediction. 




Figure 4 . Three-dimensional subsonic flow, isentropic Euler set, sharp-edged 
delta wing, 80x38x48 cell, M. * 0.3, a = 20.5°, AR = 1, e 2 = 

0 . 12 , £4 * 0.005, (a.,b.) 1 . surface pressure, 2 . static presure 
contours, Cp, 3. crossflow velocity 
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UNIFORM ROLLING IN A CONICAL FLOW 

Figure 5 shows the results for a flat plate sharp-edged delta wing which is 
undergoing uniform rolling in the counter-clockwise direction around its axis ox' at 

a constant angular speed o> = = 0.5; where w and to are the dimensional and 

L oo 

dimensionaless angular speeds, i is the wing root chord and U is the freestream 

oo 

speed. The wing angle of attack a = 0 and hence the flow is steady in the rotating 
frame of reference. Figure 5a shows the upper "l" and lower "a" surface 
pressure. As symmetric surface pressure. Figures 5b, 5c and 5d show the 
corresponding cross-flow velocity, cross-flow Mach contours an static pressure 
contours. On the cross-flow Mach contours and under the anti -symmetric leading-edge 
vortices, one notices a weak cross-flow shock. 




b. cross-flow velocity 


Fig. 5 . Uniform rolling of a delta wing, 
M =2, o=0° , (3=70°, uj=0.5 

CD 
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ROLLING OSCILLATION IN A LOCALLY-CONICAL FLOW 


The wing is given a rolling sinusoidal oscillation of u and 6 of w = - w cos kt 

e , and e = - e sin kt, where e mav = u /k and k is the dimensionless reduced 
^ nia a max o 

frequency of oscillation (K = jj— ^ - is the dimensional frequency). Choosing 

0 max = u ^ anc * % = 0*35, the corresponding k = 1.337 and the period of oscillation 

t = 4.699. Figure 6 shows the rolling oscillation motion. 



Point No. Time 9° Sense Figure 


1 

1.07 

-14.85 

CCW 

6 

2 

2.16 

-3.76 

CW 

7 

3 

3.19 

13.5, 

CW 

8 

4 

4.31 

7.46 

CCW 

9 

5 

5.35 

-11.46 

CCW 

10 

6 

6.46 

-10.63 

CW 

11 

7 

7.48 

8.17 

CW 

12 

8 

8.59 

13 

CCW 

13.1 

9 

9.69 

-5.69 

CCW 

13.2 

10 

10.7 

-14.79 

CW 

13.3 

11 

11.8 

1.03 

CW 

13.4 

12 

12.8 

14.8 

CW 

13.5 

13 

13.9 

3.93 

CCW 

13.6 

14 

15.0 

-14.01 

CCW 

13.7 


Fig. 6. Roll angle, angular speed and angular 

acceleration of the rolling oscillation 
motion. 
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ROLLING OSCILLATION (CONTINUED) 


Figure 7 shows the results for the time range t = 0 -1.07. By the end of this 

time, the wing has rolled through an angle e = -14.85° which corresponds to the end 

of the first quarter of the cycle. At t = 0, |oj| = oj 0 and |w| decreases in the 

counter-clockwise direction within that time. Figure 7a shows the surface pressure 
after each 400 time steps covering a total of 2000 time steps. On the upper surface, 
the suction pressure on the left is higher than that on the right, and the suction 
peak is moving in the spanwise direction in the positive z direction. On the lower 
surface, the surface pressure is decreasing on the left side while it is increasing 
on the right side. Figures 7b, 7c and 7d show the cross-flow velocity, cross-flow 
Mach contours and static pressure contours at t = 1.07 and e = -14.85°. At this 
moment, a large leading-edge vortex appears on the left and a small leading-edge 
vortex appears on the right. The cross-flow Mach contours show shocks above and 

below the left leading-edge vortex. It also shows the outer bow shock with varying 

strength (lower portion of the figure). 




a. surface pressure 


b. cross-flow velocity 


Fig. 7 . Rolling oscillation of a delta wing, 
M^, o=10°, P=70° , a>=0. 35, k = 1.337, 

0 max =15 ° • t=0_1 - 07 » 6“0-(-1485 0 ) 
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cross-flow Mach 


d. static pressure 


Fig. 7 . Rolling oscillation of a delta wing, 
M o> =2, a*10°, 0=70° , m- 0.35, k = 1.337, 

6 =15°, t-0-1.07, 0-O-(-14S5°) 




ROLLING OSCILLATION (CONTINUED) 

Figure 8 shows the results for the time range t = 1.07 + - 2.16. Within this 
time, the wing has reversed it direction |u|, and by the end of this time, the 
wing roll angle e = -3.76°. Figure 8a shows the surface pressure after each 400 time 
steps covering the range of time steps from 2,001-4,000. The peak suction pressure 
on the left is decreasing corresponding to a decrease in size of the left vortex 
while that on the right is increasing corresponding to an increase in size of the 
right vortex. Figures8b, 8c and 8d show the cross-flow velocity, cross-flow 
Mach contours and static-pressure contours at t = 2.16, e = -3.76°. 




b. cross-flow velocity 


Fig. 8 . Rolling oscillation of a delta wing, 

M <jo =2 , o=10° , p®70° , (0=0.35, k=1.337, 

6*=15°. t=1.07 + -2.16, e=H4.85°)-(-3.76°) 

ino a 
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ROLLING OSCILLATION (CONTINUED) 

Figure 9 shows the results for the time range t = 2.16 + - 3.19 during which the 
wing is rotating in the clockwise direction with decreasing |u|. By the end of this 
time, the roll angle is e = 13.5°. Figure 9a shows the surface pressure covering 
the range of time steps 4,001-6,000. The peak suction pressure on the left is moving 
to the left as the vortex is disappearing, and an attached flow is forming. Thepeak 
suction pressure on the right is moving inboards to the left, while the shock under 
the vortex is growing. Figures 9b, 9c and 9d show the cross-flow velocity, cross- 
flow Mach contours and static pressure contours at t = 3.19 and e = 13.5°. 


0 

-0 r«.?T7M 



a. surface pressure b. cross-flow velocity 

Fig. 9. Rolling oscillation of a delta wing, 

M oq = 2, a=10°, P=70° , to=0.35, k-1.337, 

0* =15°, t=2.16 + -3.19, 0 = (-3.76° )- 

(+13.6°) 
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c. cross-flow Mach 


d. static pressure 


Fig. 9.. Rolling oscillation of 
M =2, a=10° , 0=70°, u)=0 

(+13.5°) 


a delta wing, 
.35, k=1.337, 

0=(-3.76° )- 





ROLLING OSCILLATION (CONTINUED) 

Figure. 10(a-d) shows the results for the time range t = 3.19 + - 4.31 during 
which the winq has reversed the direction of rotation from the CW to CCW, and |to| has 
reached a zero value and then increases. The vortex on the left is growing, and the 
corresponding peak suction pressure is increasing and moving inboards to the right. 
The vortex on the right is flattening, and its peak suction pressure is decreasing 
and movino inboards to the left. By the end this time t = 4.31; the roll angle 
9 = 7.46°I 
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a. surface pressure 


b. cross-flow velocity 


Fig. 10. Rolling oscillation of a delta wing, 

14^=2, a=10° , p=70°', w=0.35, k = 1.337, 

0 max =15 °* t= 3-19 + -4.31 , 9 = (+13.5 0 )- (+7.46°) 
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ROLLING OSCILLATION (CONTINUED) 


Within the time range t = 4.31 + - 5.35, the wing is still rotating in the CCW, 
and LI has reached its maximum value and then decreases. At t = 4.7, the wing has 
already completed one cycle of oscillation. The peak suction pressure on the left is 
moving inboards to the right, and the peak suction pressure on the right is moving 
outboards to the right, Figure 11 (a--d). By the end of this time t = 5.35, the roll 
angle e = -11.46°. 
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a. surface pressure 


b. cross-flow velocity 


Fig. 11. Rolling oscillation of a delta wing, 

14^=2, a=10° , P=70°, w=0.35, k = l . 337 , 

9°° =15°, t-4.31 + -5.35, 0 = (+7.46° )- (-11.46°) 

max 


286 


ORIGINAL PAGE IS 
OF POOR QUALITY 




c. cross-flow Mach 


d. static pressure 


Fig. 11, Rolling oscillation of a delta wing, 
1^=2, a=10° , (3=70° , u=0.35, k=1.337, 

0* «15°, t=4.31 + -5.35, 0=(+7.45° )“ 

ma X 
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ROLLING OSCILLATION (CONTINUED) 

Figures 1 2a- 1 2d and 13a-13d show the results covering the ranges of time steps 
10,001-12,000 and 12 ,001-14 ,000 , respectively. Figures 5-11 cover a total of 1.56 
cycles of oscillation. The results show the successive increase, decreased and 
motion of the left and right vortices and their corresponding peak suction pressure. 
They also show the formation and disintegration of the shocks below and above the 
vortices, and the motion and strength variation of the outer bow shock. 




b. cross-flow velocity 




Fig. 12. Rolling oscillation of a delta wing, 

M oo =2, a=10° , p=70° , co=0.35 p k = l. 337, 
e” =15°, t=5.35 + -6.46, 0 = (-ll. 46°)- (-10.63°) 

ITla X 
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c. cross-flow Mach d- static pressure 


Fig. 13. Rolling oscillation of a delta wing, 

Vlj=2, a=10° , p=70°, u)=0.35, k=1.337, 

0°° =15°, t=6.46 + -7.48, 0 = (1O.63° )-(+8.17°) 
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ROLLING OSCILLATION (CONCLUDED) 

Steady state oscillation response is reached after 3 cycles of transient 
response, Figure 14 (a— g). 



Fig. 14 . Rolling oscillations from t = 7.48 
to t = 15.0, steady state oscilla- 
tion is reached. 
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Fig. 14. Continued. 
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Continued. 
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ROLLING OSCILLATION, LIFT AND ROLLING-MOMENT COEFFICIENTS 

Figure 15 shows the time history of the lift and rolling-moment coefficients 
along with roll angle variation. Steady-state oscillation response is reached by 
the third cycle, while the phase angle between C^ and ©is 90°, C M is in phase with o>. 
Although C M and e have the same frequency, C^ snows twice the value of that frequency. 
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FIGURE 15. TIME HISTORY OF THE LIFT AND ROLLING-MOMENT COEFFICIENTS 
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3. Integral Solution of Full -Potential Equation 
With and Without Embedded Euler Domains 


Formulation 


Full Potential Equation (Shock Capturing, SC; Shock Capturing-Shock Fitting, SCSF) 
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Integral Solution of Velocity Field With Explicit Shock Surface Contribution 
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IW 


Euler Equations (Integral Equation With Embedded Euler Domains) 


M a? 

St bx by e ® 

(15) 

Q • Ip, pu, pv, pej 6 

(16) 

— 2 
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( Y +l)p + (u + V )/2, h = e + p/p 
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Method of Solution 


• Shock-Capturing Shock-Fitting (SCF) Scheme: 

Shock Capturing Part 

1. In Eq. (7), Set G = q s = 0, Use Eq. (7) to Satisfy B.C.S; Eqs. (4), (6). 
Find qg and y g (piecewise linear distribution). 

2 . Calculate Initial Values of G at the Centroids of Field Elements 

(G = M* u x ). 

3. With q $ = 0, Eq. (7) is Used to Satisfy B.C.S. 

4. Calculate p and G Using Eas. (2), (3) (Type Finite Difference is Used 

for p x * p y) * 

5. Steps (3) and (4) are Repeated Until Captured Shock Location is Fixed. 
Shock Fitting Part 

6. Shock Panels are Introduced, q $ - term in Eq. (7) is activated. 

7. Eqs. (8) and (12) are Used to Find q s and p, and Eqs. (9)-(12) are Used to 

Cross the Shock. 

8. Iterative Procedure is Continued Until Convergence. 
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• Integral Equation With Embedded Euler (IEEE) Scheme (For Strong Shocks) 

1. Shock Capturing is Used to Locate the Shock. 

2. A Computational Fine-Grid Domain is Constructed Around the Shock for Euler 
Computations. 

3. With B.C. and I.C. Found From Step (1), Euler Eqs. ( 15)-( 19) are solved by Using 
a Central-Difference Finite-Volume Solver with Four-Stage Runge-Kutta Time 
Stepping and Added Second- and Fourth-Order Dissipation. 

4. Fixing q Values Found From the Euler Calculations, the Integral Ea. is Used to 
Update the B.C. 

5. The Iterative Procedure is Continued Until Convergence is Achieved. 



Fig. 16 Integral Equation Grid with an 
Embedded-Euler Domain. 
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SHOCK-FREE FLOW 


The first step to validate the computer program is to check the sensitivity of 
the IE solution to the size of the computational domain. Figure 17 shows the 
solutions for the NACA 0012 airfoil at = 0.72 and a = 0° using vortex panels only 
on the airfoil surface. We used a total°°of 140 vortex panels on the airfoil surface 
and a 64x60 field elements around the airfoil. The solutions show the surface 
pressure using two sizes of the computational domains; 2x1.5 and 3x2.5. In Figure 18, 
we repeat the same test for a lifting case of the same airfoil at M = 0.63 and 
a = 2°. The results of these two cases show that a computational domain of 2x1.5 
gives as accurate solutions as those of the 3x2.5 computational domain. 



Fig. 17 Effect of the Computational Domain Fig. 18 Effect of the Computational Domain 
Size, Surface Vortex Panels, NACA Size, surface Vortex Panels, NACA 

0012, M a =0.72, <x=0°. 0012, M =0.63, a=2° . 


300 



SHOCK-FREE FLOW 


The second numerical test is aimed at comparing the results of the IE solution 
using vortex panels only and source panels only with the solution of Euler 
equations^. Figure 19 shows the results of this test for the NACA 0012 airfoil at 
= 0.7 and a = 0°. The computational domain is 2x1.5, and the same numbers of 
surface panels and field elements as those of Figure 17 have been used. It is clear 
that the IE solution with surface vortex panels is superior to that of the source 
panels. 



o .2 .4 .6 .8 i.o 

x/c 


Fig. 19 Comparisons of IE Solution with 

Surface Vortex Panels and Surface 

Source Panels with Euler Solution, 

NACA 0012, M =0.72, ot=0° . 

00 
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TRANSONIC FLOW 


First, we present a numerical test case to show the effect of introducing the 
shock panels and their fitting as explained earlier. Figure 20 shows a comparison 
between the shock capturing results and the SCSF-scheme results for the NACA 0012 
airfoil at = 0.8 and a = 0°. It is clear that the SCSF-scheme sharpens the 

shock, as expected, with this relatively coarse grid. Next, we compare the SCSF- 
scheme with the experimental data and other computational results. Figure 21 shows 
the results of the SCSF-scheme for NACA 0012, M = 0.8 and a = 0°, along with 

comparisons with the computational results of Garabedian, Korn and Jameson 0 and the 
experimental data taken from reference 7. The SCSF-scheme took 12 iteration cycles 
of shock capturing (SC) and 13 cycles of shock fitting (SF) to achieve convergence. 
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TRANSONIC FLOW 


Figure 22 shows the results of the IFFF-scheme for the same case along with a 

comparison with the computational results of Jameson®, who also used the finite- 
volume Fuler scheme with four-stage Runge-Kutta time stepping. In the present I EEE- 
schene, the embedded Euler domain has a size of 0.5x0. 6 around the shock region with 
a grid of 25x30. This case took 10 iteration cycles of SC, 250 times cycles of Euler 
iterations to achieve a residual error of 10"'* and 5 IE cycles to update the Euler 

domain boundary conditions. 
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Fig. 22 Integral Equation with Embedded- 
Euler Domain Solution, NACA 0012, 
M =0.8, a=0° . 

CO 
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TRANSONIC FLOW 


Figures 23 and 24 show the results of the SCSF-and IEEE-schemes for NACA 64A010A, 
M = 0.796, a = 0° along with comparisons with the computational results of Edwards, 
Bland and Seidel 9 who used the TSP-equation, and the experimental data taken from 
reference 9. With the SCSF-scheme, the numbers of SC and SF iteration cycles to 
achieve convergence are the same as those of the case presented in Figure' 21. With 
the IEEE-scheme, the embedded Euler domain has a size of 0.7x0. 6 with a grid size of 
35x30. This case. Fig. 24, took 10 iteration cycles of SC, 130 time cycles of Euler 

iterations to achieve a residual error of 10”^ and 3 IE cycles to update the Euler 
domain boundary conditions. 



Fig. 23 Integral Equation Solution with Fig. 24 Integral Equation with Embedded- 

SCSF-Scheme, NACA 64A010A, Euler Domain Solution, NACA 64A010A, 

M =0.796, a=0° . M^O.796, a=0° . 
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TRANSONIC FLOW 


Figures 25 and 26 show the results of the SCSF- and IEEE-schemes for the lifting 
case of NACA 0012, = 0.75 and a = 2° along with the computational results of 

Steger and Lomax*®, and the experimental data taken from the same reference. The 
size of the grids and the number of iteration cycles used to achieve convergence are 
the same as those of the cases given in Figures 21 and 22. 




Fig. 25 Integral Equation Solution with 
SCSF-Scheme, NACA 0012, M =0.75, 
a=2° . 


Fig. 26 Integral Equation with Embedded- 
Euler Domain Solution, NACA 0012, 
M^O.75, a =2° . 
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TRANSONIC FLOW WITH STRONG SHOCKS 


For stronger shocks than those considered above the IE computational domain is 
extended in the longitudinal and lateral directions and so is the embedded Euler 
computational domain. The Euler domain is extended beyond the trailing edge to allow 
for the vorticity to be shed downstream where the overlapping region with the IE 
equation exists. Figure 27 shows a typical computational domain with details of the 
embedded Fuler domain. 



Fig. 27 Embedded Euler Domain and Grid for 
Strong Shocks. 


! 
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TRANSONIC FLOW WITH STRONG SHOCKS 


Fioure 28 shows the results of the IEEE for NACA 0012, M = 0.812 and a = 0° 

along with the experimental data of reference 7. In Figure 29, the results of the 
IEEE for NACA 0012, fl = 0.82 and a = 0° are shown along with the three-dimensional 
solution at the wing root chord of Tseng and florino , who use the IE for the TSP, 
and the results of reference 11. The size of the embedded Euler domain for these 
cases is 0.8x0. 8 and its grid size is 40x40. 
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Fig. .28 Integral Equation with Embedded- 
Euler Domain Solution, NACA 0012, 
M =0.812, a=0° . 
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TRANSOM I C FLOW WITH STRONG SHOCKS 


Figure 30 shows the results of the IEEE for NACA 0012, M m = 0.84 and a = 0° 
along with comparisons with the nonisentropic FP-solution of Whitlow, et al and 
the Euler equations solution of Jameson®. The size of the embedded Euler domain for 
this case is 1. 5x1.0 and its grid size is 60x40. This case took 10 IE iteration, 300 
time cycles of Euler iterations and 3 IE cycles to update the Euler domain boundary 
conditions. 



Fig. 30 Integral Equation with Embedded- 
Euler Domain Solution, NACA 0012, 
M =0.84, <x=0°. 

CD 
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4. Concluding Remarks 


• Two Methods have been Presented for Inviscid Transonic Flows: 

- Unsteady Euler Equations in a Rotating Frame of Reference for Transonic-Vortex 
FI ows. 

- Integral Solution of Full -Potential Equation with and without Embedded Euler 
Domains for Transonic Airfoil Flows. 

• The Computational Results Covered: 

- Steady and Unsteady Conical Vortex Flows 

- Three-Dimensional Steady Transonic Vortex Flow 

- Transonic Airfoil Flows 

• The Results are in good agreement with Other Computational Results and 
Experimental Data. 

• The Rotating Frame of Reference Solution is Potentially Efficient as Compared 
with the Space-Fixed Reference Formulation with Dynamic Gridding. 

• The Integral Equation Solution with Embedded Euler Domain is Computationally 
Efficient and as Accurate as the Euler Equations. 

• Currently the Rotating Frame of Reference Eule Solver is Applied to Three- 
Dimensional Unsteady Transonic-Vortex Flows. The IEEE-Scheme is Being Extended 
to the Unsteady Transonic Airfoil Calculations. 
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This paper reports on recent experience in calculating unsteady transonic flows by 
means of viscous-inviscid interactions with the XTRAN2L computer code (ref. 1). The 
boundary- layer method for attached flows is based upon the work of Rizzetta (ref. 2) as 
implemented in reference 3. The non-isentropic corrections of Fuglsang and Williams (ref. 
4) have also been incorporated along with the viscous interaction for some cases and initial 
results are presented. For unsteady separated flows, the inverse boundary-layer equations 
developed by Vatsa and Carter in reference 5 are used in a quasi-steady manner and 
preliminary results are presented. Currently, efforts are underway to include the viscous 
interactions in 3-D calculations in a stripwise fashion although no results for the 3-D 
work are presented herein. 


UNSTEADY TRANSONIC VISCOUS-INVISCID INTERACTIONS 


O ATTACHED FLOW (ARC, RIZZETTA) 

- THEORY WELL-DEVELOPED 

- EXTENSIVE APPLICATIONS 


O NONISENTROPIC CORRECTIONS WITH VISCOUS EFFECTS 

- INITIAL THEORY 

- LIMITED APPLICATIONS 


O SEPARATED FLOW 

-PRELIMINARY THEORY 

- LIMITED APPLICATIONS 

O 3-D FLOW 

- XTRAN3 

- CAP-TSD 
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The inviscid code used in this study is the XTRAN2L computer code described in reference 
1 . The viscous boundary layer analysis is based upon Green's lag-entrainment equations as 
described in reference 2. For attached flow, the equations are used in the direct form: 
pressure from the inviscid analysis is specified and the equations are integrated to obtain 
the boundary- layer displacement thickness 8*. For separated flows, the equations are 
inverted as described in reference 5 and the mass flow m is specified as input. In the 
inverse method, the boundary- layer displacement thickness 5* is updated using Carter's 
method (ref. 6). For both the direct and the inverse method, the effect of the viscous 
boundary layer is included in the inviscid analysis by means of the airfoil surface tangency 
boundary condition. 


BOUNDARY LAYER ANALYSIS 
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Non-isentropic modifications to the transonic small disturbance (TSD) equation were 
developed by Fuglsang and Williams in reference 4. These modifications include a 
streamwise flux that satisfies the Prandtl relations at shock jumps, convection of shock 
generated entropy in the wake, and an entropy correction in the pressure coefficient. The 
non-isentropic modifications have been incorporated into the computer code along with the 
viscous-inviscid interactions and some initial calculations are presented. 


NONISENTROPIC MODIFICATIONS TO TSD EQUATION 
(FUGLSANG AND WILLIAMS) 

O MODIFIED STREAMWISE FLUX SATISFYING PRANDTL RELATION AT SHOCK JUMPS 
O CONVECTION OF SHOCK GENERATED ENTROPY IN WAKE 
O ENTROPY CORRECTION IN PRESSURE COEFFICIENT 
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This figure illustrates one of the numerical difficulities encountered with the 
interacting boundary- layer calculation. The original version of the computer code included 
an upwind switch in the evaluation of the pressure gradient term <f>xx for input to the 
boundary- layer equations. This upwind switch introduced a discontinuity in the unsteady 
forces when the shock moved across a grid point. The dashed line in the figure shows this 
discontinuity in the moment coefficient for a typical case. The purpose of upwind switching 
in computational fluid dynamics is to account properly for the domain of dependence in the 
numerical solution of partial differential equations. However, the present application 
merely requires numerically computing the derivative of a known function. Hence, upwind 
switching is not required. The solid line in the figure shows that the moment coefficient 
varies smoothly in time when upwind switching is not used. For all results presented in 
this paper, upwind switching is not used in the boundary-layer calculation. However, the 
inviscid solution algorithm does use upwind switching in the standard manner. 


EFFECT OF UPWIND SWITCH IN 
BOUNDARY LAYER CALCULATION 

NACA 64A010A M n 0.796 
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This figure shows steady pressure distributions for the NACA 64A006 airfoil for a range 
of Mach numbers. For the subsonic cases, the viscous and inviscid results are nearly 
identical except for small differences near the trailing edge where the viscous results more 
closely match the experiment. At M = 0.850, a shock wave develops near midchord and the 
viscous calculation agrees much better with the experiment than the inviscid result. For M 
= 0.875, the shock strengthens and moves aft. In this case, both the viscous and inviscid 
calculations exhibit differences from the experiment although the viscous result is closer to 
the experiment in the vicinity of the shock. Downstream of the shock the viscous result is 
in good agreement with the experiment. 


STEADY PRESSURE DISTRIBUTION FOR 
NACA 64A006 AIRFOIL 


.8 


.6 


M = 0.800 




upper 

lower 
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ORIGINAL PAGE IS 
OF POOR QUALITY 


This figure shows plots of the boundary-layer displacement thickness for several mach 
numbers for steady calculations for the NACA 64A006 airfoil. For the two lowest values of 
Mach number, the results are subcritical and the displacement thickness increases 
smoothly in the downstream direction as Mach number is increased. At M - 0.850, the 
displacement thickness shows a slight increase due to the shock wave near midchord. For M 
= 0.875, the calculated result has a strong shock near 60% chord and the displacement 
thickness increases significantly across this shock. A further increase in Mach number to 
0.960 moves the shock off the trailing edge and the displacement thickness increases slowly 
as the trailing edge is approached. 


EFFECT OF MACH NUMBER ON DISPLACEMENT THICKNESS 
FOR NACA 64A006 AIRFOIL 



10 — 


M = 0.800 
M = 0.825 
M = 0.850 
M = 0.875 
M = 0.960 
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Unsteady pressure distributions are plotted for the NACA 64A006 airfoil with an 
oscillating flap for M = 0.850 and k - 0.242. The calculated mean pressures are similar to 
the steady pressure distributions for this Mach number. The effect of the viscous boundary 
layer can be seen in the results for the upper surface pressure distribution in the lower left 
hand side of the figure. The viscous unsteady pressure distributions agree very well with 
the experimental results, especially in the vicinity of the shock where the inviscid 
calculation shows the largest discrepancy. The interacting viscous boundary layer gives a 
substantial improvement over inviscid calculations in predicting the unsteady pressure 
distributions for this airfoil. 


UNSTEADY PRESSURE DISTRIBUTION FOR NACA 64A006 AIRFOIL 

M = 0.850 K = 0.242 




% 
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ORIGINAL PAGE IS 
OF POOR QUALITY 

This figure shows plots of the unsteady lift as a function of mach number for the NACA 
64A006 airfoil for a reduced frequency of 0.060. The results indicate that the viscous 
boundary layer corrects up to 25% of the differences between the inviscid results and the 
experiments for the lower values of Mach number. The linear theory results, also shown on 
the figure, are competitive with the CFD calculations for predicting the unsteady lift for 
most of the cases investigated for this airfoil. 


COMPARISON OF UNSTEADY LIFT FOR 
NACA 64A006 AIRFOIL 
K = 0.060 


• EXPERIMENT 

D LINEAR 

O XTRAN2L 
A PRESENT METHOD 




MACH NUMBER 
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The steady pressure distributions for the MBB-A3 supercritical airfoil at the 
supercritical design point (M = 0.765, ao = 1.5°) are plotted in this figure. The calculated 
cases are for the actual experimental values of Mach number and angle of attack rather than 
values adjusted to match flow conditions in the wind tunnel as is frequently done in 
comparison with this particular data. The experiment shows supercritical flow without a 
discernable shock wave typical of flow at the design point. The region of supercritical flow 
terminates at about x/c = 0.53. The viscous calculation indicates a moderate strength shock 
wave at nearly the same location. Away from the shock, agreement between the viscous 
calculation and the experiment is very good, although some discrepancies are noted near the 
leading edge on the lower surface. For this case, the inclusion of viscous effects yields a 
significant improvement in the calculation of the steady pressure distribution. 


STEADY PRESSURE DISTRIBUTION 
FOR MBB-A3 AIRFOIL 

M = 0.765 «0 * 1-5° 


— Viscous 

In viscid 

O Experimental upper 

□ Experimental lower 
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ORIGINAL PAGE IS 
OF POOR QUALITY 

Steady pressure distributions for the NACA 64A010A (Ames) airfoil at M = 0.796 and 
ao = 0° have been calculated with the inviscid code, the viscous interaction theory, and non- 
isentropic corrections to the viscous interaction results. As the figure shows, the viscous 
calculation agrees better with the experimental results in so far as shock location and 
strength is concerned. The non-isentropic corrections move the shock position downstream 
about 1% chord and increase the shock strength slightly. In general however, differences 
between the experiment and all three calculated results are small. 


STEADY PRESSURE DISTRIBUTION FOR NACA 64A010 AIRFOIL 
WITH NON-ISENTROPIC CORRECTIONS 
M a 0.796 


Entropy* Viscous 
Viscous 
inviscid 
C p star 

Experimental 
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This figure shows plots of the unsteady pressure distributions for the NACA 64A010A 
(Ames) airfoil oscillating in pitch at M = 0.796 and k = 0.025. In contrast to the small 
effect the non-isentropic corrections have on the steady pressure distributions, these 
corrections to the viscous interaction theory give substantially better agreement with the 
unsteady experimental results. In the vicinity of the shock wave, the modified theory 
matches the experimental points very well, whereas the inviscid calculation is quite 
different. The significant improvement of calculated unsteady pressure distributions for 
this case due to the inclusion of non-isentropic and viscous effects is particularly 
interesting because the shock wave has only moderate strength. 


UNSTEADY PRESSURE DISTRIBUTION WITH NON-ISENTROPIC 
CORRECTIONS FOR NACA 64A010 AIRFOIL 
M □ 0.796 K = 0.025 
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ORIGINAL PAGE IS 
OF POOR QUALITY 

This figure presents comparisons of unsteady lift for calculated and experimental 
results for the NACA 64A010A (Ames) airfoil oscillating in pitch. The non-isentropic and 
viscous corrections give substantial improvements in the calculated values for low values of 
reduced frequency. This is especially evident in the imaginary part of the unsteady lift. 
Neither the inviscid nor viscous calculations predict the upward trend of the experimental 
results for low values of reduced frequency. The non-isentropic and viscous calculations 
show this low frequency upward trend very well, although some discrepancies are evident 
in the mid-frequency range. This significant improvement in lift predictions for low 
frequency cases re-emphasizes the importance of non-isentropic corrections even when the 
flow field does not exhibit strong shock waves. 


COMPARISON OF UNSTEADY LIFT WITH NON-ISENTROPIC 
CORRECTIONS FOR NACA 64A010 AIRFOIL 
M = 0.796 


Experiment 

— Inviscid 

— Viscous 
Entropytviscous 
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The inverse boundary* layer code has been used to calculate several test cases in which 
the flow is separated, or close to separation, in order to demonstrate the capabilities of the 
method. This figure shows plots of the steady pressure distributions for the NLR 7301 
airfoil at M * 0.70 and ao = 3°. Note the mean angle of attack is the actual experimental 

angle of 3° and not the corrected value of 2° which is specified in the AGARD conditions for 
this case. As the figure shows, the inverse boundary -layer code predicts a pressure 
distribution which agrees reasonably well with the experiment for this very difficult case. 
The calculated shock wave is about 5% chord upstream of the experiment and slightly 
weaker. The calculation also indicates a small region of separation at the trailing edge 
whereas the experimental pressures show no evidence of trailing edge separation. 


STEADY PRESSURE DISTRIBUTION CALCULATED WITH INVERSE BL CODE 

FOR NLR 7301 AIRFOIL 

M = 0.70 ao b 3° 

— Viscous upper 

— ~ Viscous lower 

Experimental upper 
□ Experimental lower 
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This figure shows unsteady pressure distributions for the NACA 0012 airfoil oscillating 
in pitch at M * 0.599. The mean angle of attack is 4.86° and the reduced frequency is k * 
0.081. For the viscous calculation with the inverse boundary -layer method, transition is 
specified to be at 20% chord. The mean pressure distributions show the calculated results 
underestimating the suction peak near the leading edge with the inviscid calculation being 
slightly closer to the experiment in this region. Over the rest of the airfoil, both viscous 
and inviscid calculations agree well with the experiment. The unsteady pressure 
distributions on the airfoil upper surface are well predicted by both viscous and inviscid 
calculations with the viscous shock location slightly upstream of the inviscid result. The 
viscous calculation indicates that the flow is very close to separation near the maximum 
angle of attack of 7.3°. In fact, as shown in the next figure, a small change in the specified 
position of transition for the viscous calculation can result in flow separation during part of 
the oscillation cycle. 


UNSTEADY PRESSURE DISTRIBUTION CALCULATED WITH INVERSE BL CODE 

FOR NACA 0012 AIRFOIL 

M = 0.599 cto = 4.86° oi s 2.44° K □ 0.081 
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This figure shows plots of the unsteady lift and moment coefficients for the NACA 0012 
airfoil as calculated by the inverse boundary- layer method with transition specified at 10% 
chord. With transition at this location, which is upstream of the shock wave, the boundary- 
layer displacement thickness increases significantly across the shock wave. This increased 
displacement thickness causes the flow to separate just after maximum lift. The lift 
coefficient plotted in the figure clearly shows the sudden decrease in lift associated with 
flow separation and the corresponding increase in lift upon reattachment. After the flow 
reattaches, the lift coefficient continues to vary smoothly throughout the rest of the cycle. 
The inverse boundary- layer method successfully captures this flow separation and 
reattachment without difficulty. 


UNSTEADY LIFT CALCULATED WITH INVERSE BL CODE 
FOR NACA 0012 AIRFOIL 

M = 0.599 no = 4.86° ai = 2.44° K = 0.81 
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This paper has presented comparisons of calculated and experimental results for 
unsteady transonic flows over airfoils. The calculations include results obtained from a 
viscous-inviscid interaction method based upon the 2-D XTRAN2L inviscid transonic 
computer code. Non-isentropic corrections have been included in some cases. The viscous 
boundary -layer equations have been solved in the direct mode for attached flow and an 
inverse mode for separated flow. The results have demonstrated that pressures and forces 
calculated by the viscous-inviscid interaction method compare well with experimental 
results for steady and unsteady attached flows. The non-isentropic corrections with the 
viscous interaction method provided improved comparisons with experiments for unsteady 
low frequency oscillations in cases involving moderate strength shock waves. Initial 
applications of the inverse boundary- layer method have demonstrated that this method can 
calculate unsteady flow fields involving flow separation and reattachment. The results 
indicate that the viscous-inviscid interaction method can provide accurate predictions of 
viscous effects in unsteady transonic flow fields. 


CONCLUSIONS 


0 PRESSURES AND FORCES CALCULATED BY VISCOUS-INVISCID INTERACTION 
COMPARE WELL WITH EXPERIMENTS FOR STEADY AND UNSTEADY ATTACHED FLOW 


O NON-ISENTROPIC CORRECTIONS WITH VISCOUS INTERACTION YIELD IMPROVED 
COMPARISONS WITH EXPERIMENTS FOR UNSTEADY LOW FREQUENCY OSCILLATIONS 


O INVERSE BOUNDARY LAYER METHOD HAS POTENTIAL FOR CALCULATING UNSTEADY 
FLOW SEPARATION AND REATTACHMENT 
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OBJECTIVES 


This list is obvious for a symposium on transonic unsteady aerody- 
namics. The last three are the most difficult to achieve. Turbulence 
is in its usual state of affairs. 


VISCOUS-INVISCID INTERACTION 

TRANSONIC 

UNSTEADY 

TURBULENT 

THREE-DIMENSIONAL 

EFFICIENT 


ROBUST 



APPROACH 


ORIGINAL PAGE IS 
OF POOR QUALITY 


The Euler code is one we've used extensively for some time now. 
The boundary-layer code solves the three-dimensional, compressible, un- 
steady, mean flow kinetic energy integral boundary-layer equations in 
the direct mode. Inviscid-viscous coupling is handled using porosity 
boundary conditions. 


. EULER EQUATIONS 
. IMPLICIT 
. FINITE VOLUME 
. UNSTEADY 
. FLUX-VECTOR SPLIT 
. THREE-DIMENSIONAL 

. BOUNDARY-LAYER EQUATIONS 
. COMPRESSIBLE 
. THREE-DIMENSIONAL 
. UNSTEADY 

. MEAN FLOW KINETIC ENERGY 
. INTEGRAL 
. DIRECT 
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OUTLINE OF RESULTS 


This slide outlines the order of the results to follow. Steady-state 
results are considered first to validate the basic inviscid and viscous 
codes, followed by the unsteady results that have been obtained to date. 


• STEADY-STATE 

• 3-D EULER, 2-D BOUNDARY-LAYER (WING-FUSELAGE) 

• 3-D EULER, 3-D BOUNDARY-LAYER (WING) 

• UNSTEADY 

• 3-D EULER, 3-D BOUNDARY-LAYER (WING) 

• 3-D EULER, 3-D BOUNDARY-LAYER (QUASI-STEADY AIRFOIL) 

• 3-D NAVIER-STOKES (AIRFOIL) 
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ORIGINAL PAGE IS 
OF POOR QUALITY 

A two-dimensional steady version of the turbulent boundary-layer 
code was used with the three-dimensional Euler code in a strip-theory 
fashion to compute the flow about the supercritical Pathfinder wing with 
fuselage. The results are shown in Figure 1. These results are in- 
cluded simply to illustrate the type of results that might be obtained 
using the boundary-layer in a strip-theory fashion. The results were 
obtained by Dr. Keith Koenig, Mississippi State, under a NASA Langley 
grant. 
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A steady-state interactive solution using the three-dimensional 
unsteady Euler and boundary-layer codes was obtained for the ONERA M6 
wing. The streamwise momentum thickness and shape factor distributions 
at about fifty percent semi-span location are compared in Figure 2 with 
the calculations of Schmitt, Destarac, and Chavmet. An isolated experi- 
mental data point at sixty percent chord location is also shown. 


Qll^xlO 3 



Fig. 2a. Boundary layer characteristics 
on the upper surface of the 
ONERA M6 wing. M“0.84, 
Re c -11.7xl0 6 , a=3 0 . 
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Hvl -M 


Surface pressure distributions for the same solution as shown in 
Figure 2 are compared in Figure 3 with the computations of Schmitt, et 
al., and experimental data at about forty-five percent semi-span loca- 
tion. The computations of Schmitt et al. used potential flow and a 
steady state three-dimensional integral boundary-layer code. 




(a) Present calculation (b) Ref. [3] 


Fig. 3. Steady pressure distribution on the ONERA M6 Wing. M=0.84, a=3° 
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A complete three-dimensional unsteady viscous-inviscid interaction 
solution was obtained on a relatively coarse grid for the ONERA M6 wing 
as shown in Figures 4 through 6. The wing was oscillated in pitch ± 2 
degrees (mean angle of attack was 0°) about the mid-chord at a reduced 
frequency of 0.3. For this case there was little viscous effect. Fig- 
ure 4 shows unsteady viscous and inviscid surface pressure distributions 
at forty-five percent semi -span location. It is a snapshot at a - 
1 .94°. 



Fig. 4. Unsteady pressure distributions on the ONERA 

M6 wing at M“0.84, k“0.3, ct=1.94°, and y/b=0.45. 
Pitch oscillation about mid-chord -2° < a < 2°. 


ORIGINAL PAGE IS 
OF POOR QUALITY 
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ORIGINAL PAGE IS 
OF POOR QUALITY 

Figure 5 shows the phase shift of the viscous solution described on 
Figure 4. 


1.00 

0.75 

0.50 

0.Zf 

.00 

-0.Z+ 

-0.50 

-0.75 

- 1.00 


Fig. 5. Phase shift of the ONERA M6 wing lift coefficient 
(viscous solution). Pitch oscillation. 

M=0. 84, fe-0.3, -2° < a < 2°. 
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Drag Coefficient 


Figure 6 shows the viscous and lnviscid results of drag coefficient 
and number of supersonic cells for the computation described on Figure 



Period 



Fig. 6. Drag coefficient and number of supersonic cells. ONERA M6 wing pitch oscillation. 
M=0. 84 , fe-0.3, -2° < a < 2°. 
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Figure 7 is a plot of normalized lift and drag verses period for an 
NACA 0012 airfoil at M a = 0.776, Reynolds number of 23.7 x 10°, oscillat- 
ing ± 1 degree in pitch about the quarter-chord point at a reduced fre- 
quency of 0.3 using quasi-steady interaction (that is, unsteady Euler 
calculation with steady-state boundary-layer calculation). 


1-50 _ 


1.08 


0.90 


- 1.00 


- 1-50 



Fig. 7 Quasi-steady interaction for an NACA 0012 airfoil oscillating 

in pitch ± 1 degree about the quarter-chord point for M * 0.776 
and Re * 23.7 x 10° at a reduced frequency of 0.3. 
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Figure 8 is a comparison of absolute values of lift verses period 
resulting from the quasi-steady interaction solution described in Figure 
7 and an unsteady Navier-Stokes solution for the same conditions. The 
Navier-Stokes solution is courtesy of Bruce Simpson, Eglin Air Force 
Base, FL. 



Fig. 8 Lift coefficients for quasi-steady interaction and unsteady Navier- 
Stokes for an NACA 0012 airfoil oscillating in pitch ± 1 degreg 
about the quarter-chord point for M = 0.776 and Re = 23.7 x 10 
at a reduced frequency at 0.3. 


C 2 
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CONCLUSIONS 


This slide compares some of the advantages and disadvantages of 
using the Euler and boundary-layer equations for investigating unsteady 
viscous-inviscid interaction. 


ADVANTAGES 

DISADVANTAGES 

ENGINEERING ANSWERS 

. MUCH MORE DIFFICULT 

FASTER 

. SEPARATION (UNSTEADY) 

LESS STORAGE 

. DEVOTED LABOR (PARTICU- 
LAR EXPERTISE) 

GRIDDING 

. COUPLING 


. ROBUST 


. FEWER PEOPLE WORKING THE 
PROBLEM 
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INTRODUCTION 

This work describes the application of unsteady 3-D Euler and 
Navier-Stokes equations to transonic flow past rotor blades, and 
wing-alone configurations. The computer code used in this study was 
developed under the U. S. Army Research Office support at Georgia 
Tech. The transonic wing-alone calculations were supported by the 
Lockheed Georgia Company under the I RAD program. 
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OBJECTIVES 


Methods based on the transonic small disturbance theory, and the 
full potential equation have matured to a point where they may be used 
by the industries for routine aeroelastic calculations. There is now 
the need to look at higher order techniques based on the Euler and the 
Navier-Stokes equations. The higher order solvers can serve two 
purposes. First, they provide a second estimate in situations where 
potential flow theory may fail (high transonic Mach numbers, strong 
shock waves), and provide benchmark runs for the validation of 
potential flow codes. Secondly, they allow the designer to study 
phenomena such as high angle of attack transonic maneuvers, supersonic 
fighter aerodynamics, and 3-D separated flow around highly loaded 
rotor blades. 


1. To Describe a Solution Procedure for 
the Numerical Solution of the 3-D 
Compressible Viscous or Inviscid Flow 

2. Apply this procedure to a number of 
fixed and rotary wing problems of 
interest 


Figure 1 



GOVERNING EQUATIONS 


The equations governing three-dimensional unsteady compressible 
flow are the Navier-Stokes equations. If the viscous terms are 
neglected, the Euler equations result. The present solution techniques 
are designed to work efficiently with both the Navier-Stokes and the 
Euler equations. All calculations have been done on an algebraically 
generated body-fitted coordinate system (c,n»C) * which is allowed to 
move with time and follow the motion of the solid. The flow 
properties of interest at a given time level are p : the density, 
u,v,w: the velocity of the fluid in an inertial coordinate system, and 
e: the total energy of the fluid per unit volume. The quantities U,"V 
and W are the contravariant components of velocity along the 5- , n~ 
and £- directions respectively. Also, p is the pressure. 


EULER EQUATIONS 
q r + E * + F n + G $ = 0 


q 


F = J 


-1 



P 



pU 

= J- 1 

pu 



puu + £ x P 

pv 


E = J* 1 

pvu + £ y p 


p w 



pwU + £ p 


e 



(e + P)U- l t P 

pV 

puv + 7? X P 
pvV + n p 
pwV + r? 2 P 
(e + P)V-7? t P 

G = J* 1 

' pW 

PuW + $ x P 
PV W + $ y p 
pwW + $ p 
(e + P)W - l t P 


p - (7-l)(e - /2p(u 2 + v 2 + w 2 )) 


Figure 2 



GOVERNING EQUATIONS (CONTD.) 


The quantities c x » C y > 5 Z etc * are the metrics of the 
transformation computed numerically using standard second order 
accurate finite difference formulas. The quantity J is the Jacobian 
of transformation. 


NAVIER-STOKES EQUATIONS 
q r + (E-E v ) { + (F-F v ), + (G-G v )j = 0 

0 

^x r xx ^y T xy ^z T xz 

E v~ J ^x r yx "** ^y r yy ^z r yz 

£t + f t + £ t 
? x zx *y zy *z zz 

t R 4- t B + £ /3 

*x K x *y H y *z K z 



0 

Vxx + *»y T xy + Vxz 

n T + 7? T + 77 T 

'x yx *y yy 'z yz 

77 T + 77 T + 77 T 

'x zx 'y zy 'z zz 

71 R +7 ] B + 77/3 
'x^x 'y K y 'z K z 



+ £ T + £ T 
J x xx ^y xy a z xz 

K T + t f + $ T 
J x yx 5 y yy J z yz 

^x r zx ^y r zy ^z T zz 

u v + $ V 0 V + m. 



Figure 3 
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TURBULENCE MODEL 


A two-layer eddy viscosity model developed by Baldwin and Lomax 
is used in this study (Ref. 1). For the mildly separated flows 
considered here this model has proved adequate. Here 1 is the mixing 
length in the inner layer, proportional to the distance m from the wall, 
and the van Driest damping factor. In the outer layer, F max is a 
measure of the velocity scales within the shear layer, while d max 
is a measure of the length scale. At large distances from the shear 
layer the eddy viscosity is designed to approach zero through the use 
of the i ntermi ttency factor F^. 


• Baldwin-Lomax two-layer algebraic model used for 
eddy viscosity. 

• Inner Layer: 

m T = for d < d c 

= (Kd) t 1 “ exp(-d + /A + ) ] 

d + = d(pr w ) ,/2 A< 

• Outer Layer: 

= '.OlGSpc^F^ for d > d c 

F w “ min(d max F max- C 2 d mx0 2/F max > 
F max = max tV' /K ] 

F k = [1 + 5 - 5 (C3 d/d max )6rl 
Q = max|v| - min|v| 


Figure 4 
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HYBRID ALGORITHM 


In the present procedure, the time derivatives appearing in the 
governing equation are discretized using a two-point backward 
difference formula. The derivatives along the £- and £ directions have 
been kept at the new time level (n+1) where the solution is sought. 
The spanwise (n~) derivatives have been explicitly evaluated using the 
latest available information at the inboard station during odd time 
steps and outboard station during the even time step. Thus, the 
computational stencil resembles a plane Gauss-Seidel algorithm, where 
the spanwise sweeps are performed in opposite directions on alternate 
iteration levels. It may be shown from a linear stability analysis 
that this technique leads to a stable algorithm. 


Implicit Euler rule: 


r 1 - q" ♦ At|? 


n+1 


Evaluate spanwise term explicitly: 

q n+1 = q 11 - At(8|E n+1 + 6 JJ F n,n+1 + 8^G n+1 ) 

where spanwise term alternates between: 

- pH -n n+1 

i,3+l,k i r j-l,k and i,j+l,k ~ ^ijj-^k 

2 Arj 2 Zr} 


Figure 5 


TIME LINEARIZATION 


The fact that some of the quantities are to be computed at the 
new time level (n+1) means that the resulting system of equations is 
algebraic, but highly nonlinear. To avoid an iterative solution of 
non-linear equations, the well known Beam-Warming linearization (Ref. 
2) is applied to the flux terms along the and c directions. The 
result is a system of linear, block pentadiagonal equations in which 
the unknown is the 'delta 1 change in the flow properties between 
adjacent time levels. 


Second-order expansion: 

E n+1 = E n + [A n ](q n+1 - q n ) + 0(At 2 ) 

G n+1 = G n + [C n ](q n+1 - q n ) + 0(At 2 ) 

where A and C are the Jacobian matrices: 

[A] = 3E/3q and [C] = 3G/3q 

The following linear system results: 

[I + At ( 6 ^A n + 6^C n ) ]Aq Il+1 = R n ' n+1 

Aq n+1 * q 11 * 1 - q* 

R n,n+1 = - At( 5^ (E n - E J» n+1 )+6 ^ (F -F v ) n,n+1 +6 s (G n -G^ n+1 )) 


Figure 6 
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APPROXIMATE FACTORIZATION 


The direct inversion of the pentadiagonal block matrix equations 
is expensive. In literature, a number of techniques are available, 
based on the approximate factorization techniques such as L-U 
decomposition or ADI decomposition. The purpose of these techniques is 
to break up the coefficient matrix into smaller, easily inverted 
matrices. In this work an ADI factorization is used to arrive at a 
system of block tri diagonal equations, which may be inverted using the 
Thomas algorithm. 


ADI Solution in the Airfoil Plane 


Approximate factorization: 

[I + At(6 f A + 6^C) ] * [I + AtS^A] [I + At6jC] + 0(At 2 ) 

Gives two linear systems with block tridiagonal matrices: 

[I ♦ At6 { A]Aq* n+1 = R n ' n+1 
[I + At6^C]Aq n+1 = Aq* n+1 

At every time step, the solver marches through the radial 
stations explicitly performing two matrix-inversion sweeps at 
each station. 


Figure 7 
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ARTIFICIAL DISSIPATION TERMS 


The use of pure central differences to advance hyperbolic, or 
weakly parabolic equations can lead to numerical instabilities. A 
variety of artificial dissipation forms have been suggested in 
literature to overcome this instability. In this work, a fourth order 
explicit dissipation form is used, and to allow large amounts of 
explicit dissipation to be used without instability, a second order 
implicit dissipation term is added to the left side. 


• Using central differencing alone leads to odd-even 
decoupling. 

• Nonlinearities of the equations produces high- 
frequency errors which grow. 

• Second-order implicit dissipation and fourth-order 
explicit dissipation used. 

• Fourth-order implicit dissipation stabilizes the 
scheme even more, but results in penta-diagonal 
systems. 

• It can be shown that the dissipation results in an 
upwinded scheme. 


Figure 8 
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FINAL FORM OF THE DIFFERENCE EQUATIONS 

The final form of the discretized equations, used in the computer 
code are shown below. 


[I + At6gA - Ate J -1 V^A^J]Aq* n+1 = R n ' n+1 _ D fi n ' n+1 
[I + AtS^C - At ( e ( +e^) J** 1 V^A^j]Aq n+1 = Aq* n+1 

where, the explicit dissipation is: 

D e n ' n+ 1 = Ate E J _1 [(V^) 2 + (V^) 2 + (V s A s ) 2 ]jq n 

and the variable implicit coefficient is: 

+ «!><**+ 5y + £ 

As an example, 

<V £ A { ) 2 Jq" ■= (q it2 - 4q i+1 ♦ 6q A - 4^ + 

(V { A { )Jq n = (q i+1 - 2q, ♦ « i . 1 ) n 

Figure 9 
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TREATMENT OF ROTOR WAKE 


In helicopter applications, the treatment of the shed tip 
vortices requires special considerations. The finite difference grid 
is usually not large enough or fine enough to capture the many 
revolutions of the tip vortices shed by several blades. In this work, 
only the portion of the shed vorticity immediately downstream of the 
rotor blade is captured by the finite difference scheme. The rest of 
the vorticity, and the wake behind the other blades is kept track of 
using a Lagrangean approach. This approach was first proposed by Tung 
and Caradonna (Ref. 3). 


• Downwash due to tip vortex significantly affects the lift 
of the blade 

• Resolving the tip vortex by finite difference techniques 
is not possible with current computer resources 

• Effects of tip vortices lying outside of the computational 
domain must be included 



• Use CAMRAD or some other wake code to obtain 
effective partial angle of attack distribution 


Figure 10 



LIFTING VISCOUS FLOW PAST A HOVERING ROTOR 


As a first application of the solution technique described, the 
subsonic lifting flow past a two bladed rotor system in hover is 
considered. The blades were made of NACA 0012 airfoils and had a 
rectangular planform. The collective pitch was 8 degrees, and the tip 
Mach number was 0.44. There is an extensive set of experimental data 
available for this configuration (Ref. 4). Here, the computed pressure 
distrbution at a number of radial stations is plotted and compared 
with experimental data. Good agreement is observed at all locations. 



Figure 11 


363 










Lift Coefident 


SPANWISE LOADING ON A HOVERING ROTOR 

The integrated lift distributions are plotted along the rotor 
radius below. It is seen that good agreement with experiments is 
found. 

These calculations were done on a 79 x 23 x 45 grid with 50 
points at each radial location on the rotor, 11 radial locations on 
the rotor and 45 points in the normal direction. They required 3.9 
seconds per time step on the CRAY XMP. 


Lift Coefficient Distribution 



Figure 12 
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NON-LIFTING TRANSONIC FLOW OVER AN ADVANCING ROTOR 


As a second example, the transonic viscous flow past a rotor 
blade tested at ONERA is considered. The tip Mach number in this case 
was 0.6, while the advance ratio (forward speed/tip speed) was 0.45. 
In the next several pages, comparisons between experiments and the 
Navier-Stokes solutions are given at the 84% span station. For the 
sake of completeness some Euler solutions are also shown. It is seen 
that both the Euler and the Navier-Stokes solutions give acceptable 
agreement with experiments. It is also seen that the Navier-Stokes 
results predict the shock locations and strength more accurately. 
These calculations were done on a 121 x 19 x 45 grid and required 4.6 
sec per iteration on a CRAY XMP. Several thousand time steps were 
needed to advance the solution from zero degree azimuth to 360 degree 
azimuth in time, (through a full revolution). 




Figure 13 
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- 20.0 0.0 20.0 - 20.0 0.0 20.0 


LIFTING TRANSONIC FLOW PAST A RECTANGULAR WING 


As a final application of the present approach, the unsteady 
transonic flow past a rectangular supercritical wing tested at NASA 
Langley Research Center is presented. The freestream Mach number was 
0.7 and the mean angle of attack was 1.98 degrees. The wing was 
constrained to oscillate in pitch at a frequency of 10 Hz. In the 
following figures, the in-phase and the out-of-phase components of the 
surface pressure distribution are plotted at several locations. 
Overall, a reasonably good agreement is observed. 



FIGURE 14. R376E21 INVISCID REAL AMPL=1.044 2000 S/C 5000 S , M = 
0.70, a 0 = 2.0 deg. 
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- 10.0 0.0 10.0 - 10.0 0.0 10.0 


LIFTING TRANSONIC FLOW PAST A RECTANGULAR WING (concluded) 



i 1 1 1 1 1 

0.0 0.2 0.4 0.6 0.8 1.0 

X/C 


r\ = 0.310 



X/C 



X/C 


FIGURE 14 . R376E21 INVISCID IMAGINARY AMPL=1.044 2000 S/C 5000 S 

(Concluded) 
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CONCLUDING REMARKS 


A promising approach for the numerical solution of 
three-dimensional Euler and Navier-Stokes equations has been 
described. Additional work is needed for improving the efficiency of 
the present procedure. It is hoped that the techniques presented here 
will find use in fixed and rotary wing aircraft analysis. For 
additional studies and code correlations the reader is referred to 
Refs. 5-8. 


1. A solution technique for the 3-D 
Euler and Navier-Stokes equations 
has been developed. 

2. A number of interesting fixed and 
rotary wing applications have been 
presented. 

3. Additional work towards improving the 
solution efficiency is now underway. 


Figure 15 
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PURPOSE 


The purpose of the present investigation is to solve the 
Navier-Stokes equations for unsteady airfoil flows. Two primary 
types of unsteady flows are considered. The first is unsteady 
periodic flow over an airfoil at a fixed angle-of-attack past 
stall. The second is unsteady flow over an airfoil which is 
pitching either sinusoidally or with a constant-rate pitch-up 
motion. For the pitching airfoil solutions, a dynamic mesh is 
employed in the computations. All results are compared with 
experiment. 


• NAVIER-STOKES EQUATIONS APPLIED TO UNSTEADY FLOWS 
ON FIXED MESHES 

• NAVIER-STOKES EQUATIONS APPLIED TO UNSTEADY FLOWS 
ON DYNAMIC MESHES 

• INCLUDE TIME TERMS IN FLUX VECTORS 

• SINUSOIDAL PITCH SOLUTIONS 

• CONSTANT-RATE PITCH SOLUTIONS 

• COMPARISON WITH EXPERIMENT 
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GOVERNING EQUATIONS 


The governing equations used in the present analysis are the 
two-dimensional Reynolds averaged thin-layer Navi er-Stokes 
equations. They are written in generalized coordinates, with 
the n-coordinate direction along the body and the ^-coordinate 
direction normal to the body. Q represents the conserved flow 
variables. The flux vectors G and H are split according to the 
method of Van Leer , with the extension to dynamic meshes given 
by Anderson et al^. J is the coordinate transformation 
Jacobian. For an unmoving mesh, n t and are zero. 


Thin-Layer Navler-Stokes 
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FLUX SPLITTING 


Fluxes are split into a forward and a backward contribution 
according to the signs of the eigenvalues of the Jacobian 
matrices, and differenced accordingly. The split-flux 
differences are implemented as a flux balance across a cell, 
corresponding to MUSCL (Monotone Upstream-centered Schemes for 
Conservation Laws)^ type differencing. For example, the 
derivative in F at the i node in the figure can be written as 
F(Q) i + 1/2 ~ F(Q)i_ 1/2 , where eaclj F(Q) can be split into its 
forwaraand backward components F and F . State variables on 
each interface are obtained by interpolation of the conserved 
variables at the appropriate nodes. Using upwind-biasing, for 
example, conserved variables at the i-1, i, and i+1 nodes are 
used to obtain the positive contribution to F(Q) at the i+1/2 
interface. The Van Leer splitting has the advantage over more 
conventional splittings that it is continuously differentiable, 
and allows shocks to be captured with at most two interior zones. 


• Split fluxes into forward and 
backward contributions 

F(Q) = F + (Q~) + F"(Q + ) 

• Use upwind biased approximation to 
spatial derivatives 


i -1/2 i+1/2 


• 

• 

• 

i-1 

i 

i+1 


• Van Leer splitting 

• Continuously differentiable 

• Allows shocks to be captured with at most 
two (usually one) interior zones 
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NUMERICAL METHOD 


An implicit, upwi nd- bi ased , f i ni t e- vol ume scheme is used to 
numerically solve the thin-layer Navi er-Stokes equations. The 
system of equations is approximately factored and solved in two 
sweeps. The (+) and (-) superscripts indicate positive and 
negative flux split quantities. All viscous terms are centrally 
differenced. The A, B, and M matrices arise from linearizations 
of the G flux, H flux, and the viscous terms, respectively. The 
method is second order accurate in space and first order accurate 
in time. The algebraic eddy viscosity Baldwin-Lomax turbulence 
model 4 is used for all turbulent flow computations. Boundary 
conditions are applied explicitly. 


Upwind Finite-Volume Approximate-Factorization 

[ Jit * 3 n V * 3 n V ] AQ ‘ ■ - RHS 
c Jit * 3 c V + _1 a c J _1 M(c)3 aq" - 


RHS = 3 "G + + 3_ + G + 3 r H + + 3 r + H 
n n c c 

-Re _1 a^CJ -1 (t x R + £ y S) ] 

TURBULENCE MODEL: BALDWIN-LOMAX 

BOUNDARY CONDITIONS: NO-SLIP, ADIABATIC WALL ON BODY 

CHARACTERISTIC ANALYSIS IN FARFIELD 
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AIRFOIL GRID 


This figure shows a partial view of a representative grid 
used in the airfoil calculations. It is a 193 x 65 C-mesh with 
clustering in the leading and trailing edge regions. Average 
minimum normal spacing on the body is 6.4 x 10 c. The grid 
extends 30 chords from the airfoil. 


193 x 65 C-mesh 
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CONSTANT ANGLE-OF-ATTACK UNSTEADY SOLUTIONS 


The NACA 0012 airfoil was analyzed at a Reynolds number of 
3 million and a Mach number of 0.3 at several angles-of -attack up 
to and beyond stall. This figure shows computed lift coefficient 
versus angle-of -attack , alpha, in comparison With the experiments 
of Loftin and Smith” and Gregory and O'Reilly”. At 0, 10, and 15 

deg., computed values are in excellent agreement with 
experimental results. At both 16 and 18 deg. angle-of-attack , 
within the region where experiment indicates stall accompanied by 
a sudden drop in lift, the computed flowfield is unsteady and 
periodic with lift coefficients varying in the ranges indicated 
in the figure. The maximum and minimum lift values agree well 
with the corresponding experimental values before and after 
stall, respectively. The Strouhal number of the periodic flow is 
given by St = ncsina/u^ , where n is the frequency of 
oscillation, c is the airfoil chord, a is the angle-of-attack, 
and u^ is the freestream flow velocity. At 21 deg., the periodic 
oscillation is no longer present. The results indicate a nearly 
steady solution, with only a small non-periodic variation in lift 
coefficient about an average value of about 1.05. 


NACA 0012 airfoil, Re - 3 million 
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STEADY VORTICITY CONTOURS 


Computed vorticity contours are shown for the NACA 0012 
airfoil at 15 and 21 deg. angle-of -attack . At 15 deg., prior to 
stall onset, the vorticity is concentrated in a relatively thin- 
layer near the airfoil surface and behind the trailing edge. At 
21 deg. the airfoil shows a region of massive separation above 
the airfoil upper surface. Vorticity levels are much higher than 
the 15 deg. case, with the strongest clockwise vorticity 
concentrated near the leading edge on the upper surface and the 
strongest counterclockwise vorticity again behind the trailing 
edge . 


NACA 0012 airfoil, Re - 3 million 



UNSTEADY VORTICITY CONTOURS 


At 16 deg. angle-of-attack , the flowfield oscillates in a 
periodic manner, with the lift varying between a minimum of 0.89 
and a maximum of 1.60. Vorticity contours are shown at four 
points in the unsteady periodic cycle. The cyclic nature of the 
flowfield is characterized by the shedding of a leading edge 
vortex near maximum lift. 


NACA 0012 airfoil, Re > 3 million, a - 16° 
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FORCED PITCHING SOLUTIONS 


Computations have been performed in which the airfoil 
undergoes a forced pitching motion about its quarter chord. Two 
types of motion have been explored and results compared with 
experiment. The first is a sinusoidal pitching motion at high 
Reynolds number (turbulent flow) and the second is a constant- 
rate pitch from 0 to 60 deg. angle-of-attack at Re = 45,000 
( laminar flow) . 


a( x) 
a( x) 


k = reduced frequency * uc/u^ 
a) = frequency (rad/sec) 
c = chord 

u w = freestream velocity 
t *° = time, nondimensionalized by c/a o 
a = freestream speed of sound 

oo 


= + a- sin(M kx) sinusoidal pitch 

0 1 00 

= a rt + M kx constant-rate pitch 

o ®° 
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SINUSOIDAL PITCH SOLUTION WITH NO STALL 


The NACA 0012 airfoil was analyzed at Re = 3.89 x 10 , 

M = 0.301, a- = 7.97 deg., a 1 = 4.91 deg., and a reduced 
frequency k = 0.398. This corresponds with Case 7111 from 
McCroskey et al^ . Transition to turbulence is fixed at the 
leading edge in the computations but was not fixed in the 
experiment. The time step taken for the computations is 0.05. 
This figure shows the lift and moment coefficients as a function 
of angle-of-attack . The thin-layer Navier-Stokes code predicts a 
shallower lift versus alpha slope, slightly o ver pr ed i ct i ng the 
minimum lift and underpredicting the maximum lift. Lift values 
for increasing alpha are on the upper half of the lift curve for 
both theory and experiment. The moment coefficient is in good 
agreement with experiment when angle-of-attack is increasing 
(lower half of the curve), but underpredicts the moment when 
angle-of-attack is decreasing. 


NACA 0012 airfoil, Re - 3.89 x 10 6 , 

M - 0.301, a Q - 7.97°, o 1 - 4.91°, k - 0.398 
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PRESSURE COEFFICIENTS 


Pressure coefficients are compared with experiment at six 
times in the pitching cycle for Case 7111. Overall, the present 
method does well in predicting the shapes of the pressure curves 


M 


NACA 0012 airfoil, Re - 3.89 x 10 6 , 
0.301, o Q - 7.97°, <*1 - 4.91°, k - 0.398 
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SINUSOIDAL PITCH SOLUTION WITH DEEP STALL 


A solution with deep stall was computed at Re = 3-76 x 10^, 
M = 0.292, a Q = 14.84 deg., a 1 - 9.88 deg., k = 0.202. These 
conditions correspond with Case 14210 from McCroskey et al?. in 
this case, both theory and experiment fixed transition to 
turbulence at the leading edge. Time step for the computations 
is 0.05. As seen from plots of lift and moment coefficient 
versus angle-of -attack , the computations agree well with 
experiment only as the angle-of-attack is increasing (upper 
portion of the lift curve, lower portion of the moment curve). 
However, theory predicts stall later than experiment. At all 
other points in the pitching cycle, theory and experiment are 
only qualitatively similar. Theory shows oscillations, which 
are not present in the experiment, in the lift and moment curves 
as the angle-of-attack decreases. 


M 


NACA 0012 airfoil, Re 
0.292, ccq - 14.84°, 


3.76 x 1 0 6 , 

9.88° , k - 0.202 


4T 


3r 


o~2b 


1b 


oL 

0 


Thin-layer N-S 

D Experiment 



Thin-layer N-S 

D Experiment 



388 


ORIGINAL page is 
OF POOR QUALITY 

AGARD CASE SINUSOIDAL PITCH SOLUTIONS 


Computations were performed to compare with two unsteady 
solutions from Landon . The first, corresponding with AGARD 
Case 2, is at Re = 4.8 x 10 , M - 0.6, a Q = 3.16 deg., 
a. = 4.59 deg., k = 0.1622. The second corresponds with AGARD 
Case 3: Re = 4.8 x 10 , M = 0.6, = 4.86 deg., « 2.44 deg., 

k = 0.1620. Lift and moment coefficients for each case are shown 
in the figure. In both cases, for theory and experiment, the 
lift values for increasing alpha are on the lower portion of the 
curve. The computations do fairly well to predict lift at the 
lower angles-of-attack but underpredict the lift. at the high end 
of the cycles. Moment coefficients are under predi cted everywhere 
in the cycles, although for alpha increasing (lower portion of 
the curve) results are in closer agreement than for alpha 
decreasing. Results for Case 3 are very similar to those 
obtained by Howlett 9 using a small-disturbance potential code 
coupled with an inverse boundary-layer method. Howlett found 
that results near maximum lift are highly sensitive to the 
transition location. His results shown here have transition set 
at 20$. The thin-layer Navier-Stokes computations set transition 
at the leading edge, while experimental transition is 
unspecified. 


AGARD Case 2 


AGARD Case 3 


NACA 0012 airfoil, Re - <1.8 x 10 6 , M - 0.6 




a Q - 3.16*. a, - 4.59*. k - 0.1622 


a 0 
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FORCE COMPONENTS AND PHASE LAG 


The real and imaginary components of lift and drag were 
computed for AGARD Cases 2 and 3 and compared with experiment. 
Agreement is fair, with the largest discrepancies in the 
imaginary component of moment for both cases. The 
values a. and a are the phase angles by which lift and moment 
cycles lead or Tag the angle-of -attack cycle. 
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CONSTANT-RATE PITCH SOLUTIONS, u> = 460 DEG/SEC 


At a Reynolds number of 45,000, laminar flow computations of 
an NACA 0015 airfoil pitched up at a constant rate of 460 deg. /sec 
(k « 0.2007) are compared with smoke wire flow visualizations of 
Helin and Walker 10 . Three angles-of-attack are shown in the 
figure. Flow is from right to left. With a time step of 0.05, 
computed velocity vectors show the same general trend as 
experiment, although the center of the shed leading-edge vortex 
does not convect downstream as quickly in the computations as it 
does in the experiment. At both 45 and 60 deg., there is 
reversed flow over most of the airfoil upper surface, due to the 
shed vortex. 


NACA 0015 airfoil, Re » 45,000, laminar flow 

« - 460 deg . /sec 



EXPERIMENT 
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CONSTANT-RATE PITCH SOLUTIONS, w 


1380 DEG/SEC 


At a higher rate of pitch (k - 0.6021), computation with a 
time step of 0.02 shows a leading-edge vortex growth rate in good 
agreement with experiment. This vortex is much smaller in size 
than that for the lower pitch rate of 460 deg. /sec. However, the 
computational analysis does not show the second region of 
separated flow near the trailing edge that is seen in the flow 
visual i zati ons . 


NACA 0015 airfoil, Re - 45,000, laminar flow 
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CONCLUSIONS 


A finite volume implicit approximate factorization method 
which solves the thin-layer Navier-Stokes equations has been used 
to predict unsteady turbulent-flow airfoil behavior. At a 
constant angle-of -attack of 16 deg., the NACA 0012 airfoil 
exhibits an unsteady periodic flowfield with the lift coefficient 
oscillating between 0.89 and 1.60. The Strouhal number is 
0.028. Results are similar at 18 deg., with a Strouhal number of 
0.033. A leading-edge vortex is shed periodically near maximum 
lift . 

Dynamic mesh solutions for unstalled airfoil flows show 
general agreement with experimental pressure coefficients. 
However, moment coefficients and the maximum lift value are 
underpredicted. The deep stall case shows some agreement with 
experiment for increasing angle-of -attack , but is only 
qualitatively comparable past stall and for decreasing angle-of- 
attack. Laminar-flow computations of a constant-rate pitch-up 
NACA 0015 airfoil show that increasing pitch rate slows 
separation. Computed velocity vectors agree qualitatively with 
experimental flow visualizations. 


• CONSTANT ANGLE-OF-ATTACK 

• NACA 0012 FLOWFIELD UNSTEADY AT a - 1 6 
AND 18 DEG. 

• lift oscillates within range of 
experiment 

• St = 0.03 

§ PERIODICALLY SHED LEADING-EDGE VORTEX 

• SINUSOIDAL PITCH 

• UNSTALLED CASES 

• pressure coefficients agree in 

general with experiment 

• moment coefficients and maximum 
lift underpredicted 

• DEEP STALL AGREES ONLY FOR a INCREASING 

• CONSTANT-RATE PITCH 

• . J INCREASING PITCH RATE LESSENS SEPARATION 

• • QUALITATIVE AGREEMENT WITH EXPERIMENT 
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OBLIQUE-WING RESEARCH AIRCRAFT PROGRAM 

The advantages of oblique-wings have been the subject of 
numerous theoretical studies, wind tunnel tests, low speed flight 
models, and finally a low speed manned demonstrator, the AD-1 (ref. 

1) . An oblique-wing configuration is well suited for a Navy fleet 
defense mission and a supersonic transport (Mach <1.6). An 
excellent review of the historical development of oblique— wing 
technology is presented in reference 2; references 3 and 4 discuss 
potential applications. NASA's Oblique-Wing Research Aircraft 
(OWRA) program is directed at the development and flight test of a 
full scale supersonic demonstrator which will address the key 
technological challenges. The specific objectives of the OWRA 
program are 1) establish the necessary technology base required to 
translate theoretical and experimental results into practical, 
mission oriented designs, 2) design, fabrication and flight test an 
oblique-wing aircraft throughout a realistic flight envelope, and 
3) develop and validate design and analysis tools for asymmetric 
aircraft configurations. 



Objectives 

• Establish a technology base 
for oblique wing concepts 
which can be applied to 
mission-oriented aircraft 
designs 


• Design, fabricate, and flight 

test an oblique wing through- 
out a realistic flight envelope 

• Develop and evaluate design 

tools for asymmetric aero- 
dynamic configurations 
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OBLIQUE -WING AERODYNAMIC ADVANTAGES 

Theoretical aerodynamic advantages of oblique wings have been the 
subject of numerous studies over the years . The variable sweep 
aspect of course permits optimization with Mach number thus 
yielding efficient flight for the subsonic cruise/loiter condition 
while also providing for efficient supersonic dash/cruise 
capability. As shown in the figure the induced drag is minimized 
for a zero sweep, maximum aspect ratio condition; this advantage is 
independent of symmetrically swept or obliquely swept aircraft . In 
the supersonic regime, the oblique type wing has a significant 
advantage (over a symmetrically swept wing) in that it produces 
less wave drag since the wing volume is distributed over a greater 
length . 
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OBLIQUE-WING AERODYNAMIC CENTER SHIFT 


An oblique-wing configuration also provides a major advantage in 
that sweep does not produce an aerodynamic center shift as does a 
symmetric swept configuration. This minimizes trim drag penalties 
due to aerodynamic center shifts and reduces tail loads, thus 
resulting in a lighter structure and eliminating center-of-gravity 
control as a function of wing sweep. 
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OBLIQUE-WING STRUCTURAL ADVANTAGE 


An oblique-wing configuration has a number of significant 
structural advantages over a symmetrically swept wing configuration, 
the most obvious being a single pivot requirement. A single pivot 
results in both cost and weight savings and other factors that 
accrue from maintaining one as opposed to two pivots. On an 
oblique-wing, the lift forces pass essentially through the center 
of the pivot independent of sweep angle, thus minimizing bending 
and torque loads transmitted through the pivot. For symmetrical 
swept configurations, offset lift forces produce significant 
bending and torque forces transmitted through the pivot, which in 
turn requires a 'beefed-up' pivot/substructure assembly and results 
in a major weight penalty. 


Symmetric wing sweep 



• Single pivot 
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OWRA UNSTEADY PRESSURE MEASUREMENT EXPERIMENT 


Efficient aircraft design increasingly relies on predictions; 
therefore, in an attempt to improve vehicle aeroelastic design and 
prediction techniques, an experiment will be implemented on the 
OWRA which will measure unsteady pressures. The unsteady pressure 
survey will use remote sensing (pneumatic lines) to measure 
pressures on thirteen chords covering the full span of the wing; 
each chord will consist of approximately 30 orifices and will be 
sampled 400 times per second. The approach is similar to that 
reported for the laboratory experiment of reference 5 and the wind 
tunnel test of reference 6. A limited number of in situ 
measurements will be taken and used to correct the pneumatic 
measurements for magnitude and phase. Controlled data will be 
gathered using preprogrammed aileron excitation algorithms. The 
data base will be used for correlation with currently used unsteady 
aero codes and will also provide a valuable data base for 
evaluation of future codes. It is anticipated that the unsteady 
pressure measurements will prove valuable in analysis of other 
unique flow phenomena and provide insight into effects such as 
vortex flow patterns and vortex and/or shock induced oscillations 
should they occur. 


\ 


Objectives: 

- Develop unsteady pressure data base 

Full span 

13 chords; 30 orifices/chord 
400 samples per second 

- Correlation with current unsteady codes 

- Data base for future code development 
~ Identify unique flow phenomena 

Vortex flow 

Vortex / shock induced oscillations 
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OWRA UNSTEADY PRESSURE SENSOR LAYOUT 


A pressure sensing system will be implemented on the OWRA which 
will be capable of measuring both static and unsteady pressures. 

The primary system will acquire data pneumatically using remotely 
located electronically scanned pressure (ESP) modules located 
either in front of the forward spar or aft of the rear spar. This 
arrangement will provide a cost effective and readily maintainable 
system since access to the wing box will not be possible after it 
is sealed. The pneumatic system will consist of equal lengths of 
tubing connecting the orifices and the transducer. Current plans 
call for use of approximately 4 foot lengths of 0.060 " ID tubing. 
The data will be corrected based on in situ unsteady pressure 
measurements made at two chord locations and one in situ 
measurement made for each of the other chord locations . At the 
maximum skew angle of 65 degrees, the left wing overlays the left 
horizontal stabilizer and as such, leads to interesting aerodynamic 
interactions. In order to assist in analysis of this effect, 
unsteady pressures will also be measured for two horizontal tail 
chord locations. 
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OWRA MODAL RESPONSE SURVEY LAYOUT 


Correlation of predicted and experimental unsteady aerodynamics 
requires an accurate mode shape of the wing. The figure below 
illustrates the planview layout for accelerometers used for 
defining the wing mode shape. In addition, these accelerometers 
will in general meet the requirements for flutter and 
aeroservoelastic stability clearance work. The unsymmetric nature 
of the OWRA leads to unsymmetrical leading edge suction forces 
which could in turn develop significant in-plane wing motion. 
Therefore triaxial accelerometers will be located along the wing. 
Additional accelerometers will also be located on the fuselage to 
identify later bending and torsional characteristics. The sample 
rate of the accelerometers will be identical to the rate used for 
unsteady pressure measurements. 
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F-15 UNSTEADY PRESSURE EXPERIMENT 


A validation of the unsteady pressure measurement system proposed 
for the OWRA was conducted on an F-15 experimental aircraft at 
Dryden . The validation consisted primarily of a parametric 
evaluation of line length and orifice/tubing diameters. An 
auxiliary objective was to demonstrate that the ESP module could be 
driven and data recorded at 500 samples/sec in a flight 
environment with no adverse effects on data quality. The 
experiment consisted of orifices located at 10% chord, one-half 
inch apart, and approximately mid-span on the upper surface of the 
right wing of the F-15. The orifice/tubing inside diameters 
evaluated were 0.020, 0.040, and 0.060 inches with tubing lengths 
of two, four, and eight feet being changeable between flights. An 
in situ measurement consisted of a 0.060 inch orifice/tubing 
diameter connected to the same ESP module but with a minimal line 
length, six inches. Flights to date have obtained excellent 
quality data for both two and four foot line lengths. 


V 
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F— 15 UNSTEADY PRESSURE SYSTEM SCHEMATIC 


An end-to-end schematic of the instrumentation system used in the 
F-15 experiment is presented in the figure below. The ESP module 
consists of 32 flush-diaphragm, strain-gage-type 

differential-pressure transducers. The reference side of the ESP 
module is connected to an ambient pressure reservoir which is 
vented to the interior wing cavity. The purpose of the reservoir 
is to allow the reference pressure (backside of the ESP 
differential transducer) to adjust to changes in altitude without 
any high frequency pressure oscillations. The module multiplexes 
the individual port measurements with the output routed to a 10 bit 
PCM system and recorded on on-board tape. There is no signal 
conditioning of the individual ESP port data prior to 
multiplexing. The ESP module transducers are ranged for +/- 5 psi. 

The ESP module is operated and data was recorded at 250 sps for 
the first few flights and subsequently increased to 500 sps. A 
heater blanket was installed to maintain a constant temperature on 
the ESP module. A check on the quality of the PSI transducer is 
obtained by plumbing one of the ports directly to the reservoir. 

The absolute transducer on the reference can allows for the 
absolute chordwise pressure measurements as well. A factor which 
makes the ESP transducer outstanding for dynamic pressure 
measurements is its minimal internal volume. The internal diameter 
of each ESP port is 0.040 inch with no increase in diameter at the 
diaphragm face and as such the transducer volume can be 
analytically modeled as a 0.040 inch ID line length extension to 
the orifice connecting tubing. 
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F— 15 UNSTEADY PRESSURE TRANSFER FUNCTION ANALYSIS 


A transfer function analysis was performed on approximately one 
minute of data obtained in a moderate g wind up turn at a high 
subsonic Mach number. Using the 0.060" tubing as the reference, both 
the 0.040" and 0.020" tubing show little attenuation to at least 10 
Hz although the 0.020" tubing does attenuate at a much more rapid 
rate than the 0.040" tubing once the break frequency is past. The 
poorer characteristics of the 0.020" tubing are also indicated by 
its significantly worse phase angle. Although not plotted on this 
figure, the coherence of the two transfer functions was also 
determined. For the 0.040" transfer function, the coherence starts 
at one (perfect correlation) for low frequencies and gradually 
decreases to one-half (reasonably good correlation) at 100 Hz. The 
0.020" transfer function also starts at a coherence of one for low 
frequencies but degrades to a value of zero (no correlation) at 
approximately 80 Hz. 
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OWRA FLUTTER MODEL TEST 


Unique aerodynamic characteristics of oblique-wing configurations 
have the potential for producing unusual flutter type 
characteristics and/or other instabilities. A wind tunnel flutter 
model test will be performed in the TDT to both provide data for 
validation of aeroelastic analysis codes prior to first flight and 
to support an efficient and rapid envelope clearing process. In 
order to maximize the return on the test, the model will be 
designed to flutter (or encounter some other type of instability 
unique to oblique-wing configurations) within the tunnel. 
Identification of transonic flutter/instability characteristics is 
of prime importance. There are tentative plans to obtain limited 
unsteady pressure measurements for both code validation and 
correlation with flight results. Preliminary studies have been 
performed to identify critical DOF for flutter model tests of 
oblique configurations. An 'oblique' mode has been identified with 
a 5 DOF model which still retains its characteristics with the 
three rotational DOF's. 


Problem: Unique, unsymmetrical configuration presents 

potential for unusual flutter type instabilities. No 
flight experience is available for oblique-wing 
configurations. 

Objectives: 

- Design model to flutter; test in TDT 

- Identify transonic flutter characteristics 

- Correlate with predictions 

- Identification of unique instability phenomena 

- Limited unsteady pressure measurements 

Status: 

- Identify important DOF 

- Preliminary study has identified an 'oblique' mode 
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FLUTTER/AEROSERVOELASTIC CODE VALIDATION 


An interdisciplinary analysis code (STARS), which is capable of 
performing flutter and aeroservoelast ic analyses, has been 
developed. The structures module has a large library of elements 
and in conjunction with numerical analysis routines, is capable of 
efficiently performing statics, vibration, buckling, and dynamic 
response analysis of structures. In order to accommodate 
unsymmetr ical supersonic conditions, the potential gradient method 
(PGM) unsteady aero code of Appa is being implemented into the aero 
module of STARS; subsonic unsteady aero code will continue to be 
doublet lattice. Linear flutter models are developed and transformed 
to the body axis coordinate system and are subsequently augmented with 
the control law. Stability analysis is performed using hybrid 
techniques. The major research benefit of the OWRA program will be 
validation of design and analysis tools. As such, the structural model 
will be validated and updated based on ground vibration test (GVT) 
results. The unsteady aero codes will be correlated with 
experimentally measured unsteady pressures. 


STARS: In-house analytical code 

- Specialized structural modeling 

- Efficient matrix manipulation 

- Implement PGM code 
Validate structural model 

- Fuselage GVT 

- Wing GVT 

- Complete A/C GVT 

Validate unsteady aero code with flight data 
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HIGH ANGLE OF ATTACK AERODYNAMICS 


As angle of attack increases, the F-8 OWRA will exhibit 
non-linearities in all flight axes. At high wing sweeps the 
increase in spanwise flow and the formation of a leading edge 
vortex can occur at relatively low angles of attack (6 to 8 
deg) . Because of the asymmetry of the vehicle these effects 
will not be balanced in the lateral directional axis. At 
higher angles of attack, regions of spanwise flow also form in 
an asymmetric pattern, generally progressing from the trailing 
wing tip. In addition to these characteristics, which will 
effect the vehicle flight dynamics, other unusual features have 
been observed in water tunnel studies such as the interaction of 
parallel spanwise vortices on the leading wing panel. Further 
water tunnel studies will be conducted this summer to document 
the flight configuration and to note distinctions between the 
various preliminary design planforms. A comprehensive 
computational fluid dynamics (CFD) analysis is underway at 
Ames-Moffett to develop a Reynolds averaged Navier-Stokes 
solution of the complete vehicle. Preliminary results of 
the wing alone at an angle of attack of 10 deg show good 
correlation of the spanwise flow and vortex formation with 
the water tunnel results. During the flight program, unsteady 
pressure data will be used to identify vortex flow and regions 
of separated flow. The vehicle is also equipped with a tail 
mounted camera which can be used for tuft studies and smoke 
flow visualization. Flight measurements of the vehicle 
forces and moments will be used for correlation with the 
flow visualization results. A similar correlation was 
made during the AD-1 flight program. 


Characteristics 

Significant spanwise flow 
Strong spanwise vortices 
Asymmetric regions of separated flow 
Potential dynamic interaction of vortices 

Data Sources 

Water tunnel studies 

Reynolds averaged Navier Stokes solutions 
Flight testing 
Unsteady pressures 
Flow visualization, tufts and smoke 
Measured vehicle forces and moments 
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OWRA FLIGHT DEFLECTION EXPERIMENT 

In order to support research activities on the OWRA, an accurate 
determination of the wings deflected shape in flight is required to 
validate the wing stiffness and load distribution predictions which, 
because of the wings unconventional attitude, could produce some 
unpredicted pressure distributions. The deflections will also be 
used for definition of in-flight shape for correlation of CFD codes 
with flight determined static pressure distributions. The 
electro-optical system to be used has been developed at Dryden and 
used quite successfully on both the HiMAT and X-29 aircraft. 


OBJECTIVE 

Evaluate the ability of analytical codes to predict structural 
loads and deflections and pressure distributions 


APPROACH 

Measure in-flight deflections to correlate with predictions 

Define in-flight shape for correlation of pressure data 
with CFD codes 
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OWRA FLIGHT DEFLECTION MEASUREMENT SYSTEM 

The flight deflection measurement system is electro-optical and 
consists of wing located targets emitting a light source which is 
in turn received by a wing center line located receiver. The 
system to be implemented on the OWRA will have the targets located 
on the wing upper surface just over the wings fore and aft spars . 
The system is capable of measuring both bending deflections and 
twist angles . Sixteen targets per wing semi-span will be 
installed; the targets have a spherical section shape with a base 
diameter of approximately 1.5 inches and a height of 1 inch. The 
receiver (located at the wing pivot) will be housed in an 
aerodynamically shaped blister with a height of approximately six 
inches and slightly larger base dimensions. The entire deflection 
measurement system will be removable so as not to interfere with 
either the static or unsteady pressure experimental data. 






RECEIVERS 

TARGETS 
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ELECTRO-OPTICAL FLIGHT DEFLECTION MEASUREMENT SYSTEM 


The major elements of the deflection measurement system are a 
control unit, two receivers, two target drivers, and the targets. 

The targets house an LED which is turned on and off sequentially 
as a command is cycled through the various targets. Light from the 
target LED is then sensed by the receiver and is focused on a light 
sensitive diode array. The signal produced by the diode array is 
proportional to the wing deflection. This signal is sent to the 
control unit which, in turn, sends it to the PCMQ system for 
recording on magnetic tape or telemetering to a ground station for 
real time display. The control unit contains all of the measuring 
logic for operation of the system, and the target driver serves as a 
relay in providing from 2 to 5 amperes of pulsed current to each 
target. The system has a resolution of approximately 0.03 inch 
at a 10 foot range and a sample rate of approximately 7 sps for 
each target . 


System Diagram 
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FLIGHT DETERMINED DEFLECTION MEASUREMENT RESULTS 

The electro-optical deflection measurement system has been used 
successfully to obtain flight results on both the HiMAT and X-29A 
research aircraft as illustrated in this figure. The X-29A flight 
measured wing panel streamwise twist distributions are derived from 
front and rear spar deflections measured during a 5 g wind up turn 
maneuver and are compared to calculated twist data. Also shown are 
HiMAT flight measurements compared to NASTRAN calculated data for 
the 8 g maneuver design point. The deflection measurements of the 
HiMAT wing played a major role in evaluating the performance of the 
aeroelastically tailored composite wing. 


Flight Measurements 

X-29A 5-g maneuver HiMAT 8-g maneuver 
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OWRA SCHEDULE 


The preliminary design phase of the project is complete and has 
resulted in a wing configuration for which construction is ready to 
be initiated. The wing area of the current configuration is 300 
sq. ft. and is a 50% increase over the feasibility study design. 
Wing final design and construction will begin this summer and be 
complete in mid 1990. Testing of the flight configuration force 
and moment model will begin this fall and upon completion, work can 
be initiated on the rigid body flight control system. Upon receipt 
of the wing, extensive ground testing will be conducted to verify 
loads and dynamic structural modeling. Extensive systems checkout 
will be conducted with emphasis on the wing and fault tolerant 
processor (FTP) based flight control system. A first flight is 
anticipated for early 1991. 


Wing construction - July 1987 - July 1990 
Systems checkout - July 1990 - April 1991 
First flight - May 1991 
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INTRODUCTION 


One aircraft configuration that shows great promise in 
achieving high performance is that of an asymmetrically swept 
wing, shown in Figure 1 (ref. 1, 2). When compared to 
conventional swept wings, these advantages include higher lift-to- 
drag ratios and reduced takeoff and landing speeds, which 
translate into increased performance in terms of fuel consumption, 
loiter time, range, etc. However, the oblique wing has a number 
of disadvantages because of its asymmetric configuration. 
Referring to Figure 1, consider the swept oblique wing shown to 
have an upward bending deflection, such that lines AB and A'B' 
represent lines of constant upward bending displacement. For the 
aft-swept portion of the wing, the airflow will see line CB. 
Since point B deflects upward more than point C (due to the 
bending displacement increasing from the wing pivot to the wing 
tip), the airflow will see a downward twist along CB due to the 
bending displacement. This bend-up/twist-down phenomenon is 
referred to as "wash-out”. The forward-swept wing, on the other 
hand, will have the airflow seeing a nose-up twist due to bending 
since point C' deflects more than point B'. This bend-up/twist-up 
is called "wash-in". The increase in angle of attack associated 
with wash-in will increase the wing load, which will tend to 
increase the bending deflection and hence wash-in twist even 
further. Thus, divergence becomes a concern with the forward- 
swept wing (e.g., the X-29). Also, because the two portions of 
the wing undergo different bend/twist behaviour, the swept oblique 
wing will have a roll imbalance due to the different loadings on 
the forward- and aft-swept portions of the wing. 

The question is, then, how to best achieve maximum stability 
and roll equilibrium without compromising performance. Using 
aeroelastic tailoring to enhance aeroelastic stability and control 
has been demonstrated in several analyses, especially for the 
forward-swept wing (ref. 3, 4, 5). Since the oblique wing has a 
forward-swept half, aeroelastic tailoring is also potentially 
beneficial for an oblique wing design. For a basic discussion of 
aeroelastic tailoring, see references 6 and 7. 
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STATIC AEROELASTIC COMPUTATIONAL PROCEDURE 


The static aeroelastic computational procedure (Figure Z> 
was developed to study the basic effects of aeroelast ically 
tailoring an oblique wing through the use of composite materials. 
First, the geometry is defined for the oblique wing, which may 
have deflected control surfaces. In this analysis the oblique 
wing model has two outboard ailerons deflected an equal but 
opposite amount as input by the user of the computational 
procedure. This geometry is then submitted to the full potential 
code FL022 for aerodynamic analysis (ref. 8, 9, 10). The output 
is a pressure distribution over the wing. After the pressure load 
has been converted to equivalent loads P, an equivalent plate 
program, developed by Dr. Gary Giles at NASA-Langley , is invoked 
for each half of the oblique wing (ref. 11). From the 
structural definition of the wing (input by the user) and the 
equivalent loads, the plate program calculates a set of 
coefficients C, from which the displacement of the wing due to the 
aerodynamic loads is defined in polynomial form. The wing shape 
is then deflected according to the calculated displacement. This 
deformed wing geometry is then input to FL022 for aerodynamic 
analysis, and the aeroelastic procedure is repeated until a 
converged deformed shape has been obtained for the flexible 
composite wing. Usually only 3-4 aeroelastic iterations are 
required before a converged shape, i.e., a shape consistent with 
the aerodynamic loads calculated by FL022, is reached. 
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Figure 2 
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APPLICATION OF AEROELASTIC TAILORING 


Each aeroelastic computational run described above involves 
one aileron deflection and one value of the amount of aeroelastic 
tailoring applied to the wing. This application of aeroelastic 
tailoring is achieved in this analysis by simply rotating the 
basic composite skin laminate of the wing by an angle 0 (see 
figure below). Recall that a swept oblique wing exhibits a roll 
imbalance. If asymmetric composite tailoring is applied to the 
wing, i.e., the aft-swept half of the wing is given a wash-in 
structure to counteract its wash-out twist due to bending (recall 
figure 1), and the forward-swept half is given a wash-out 
structure to alleviate its wash-in twist due to bending, the 
oblique wing will aeroelastically desweep in that it will 
aeroelastically behave as if the wing had less sweep. This is 
desired since an unswept oblique wing does not have a roll 
imbalance. Thus, asymmetric tailoring could alleviate the roll 
problem of the oblique wing by an aeroelastic desweeping, while 
the wing would still retain the aerodynamic advantages of being 
swept. The tailoring is simply applied by rotating the composite 
laminate an angle 6 as shown below. (Figure 3.). The wing can thus 
be trimmed in roll with aileron deflection or asymmetric tailoring , 
or a combination of both, as seen next. 


APPLICATION OF WASH-OUT/WASH-IN TO COMPOSITE WING 


Air flow 


Dominant fiber direction 
of skin laminate (90°) 



Forward wing 
[90°/+45°/0°] Laminate 


Forward wing 

[[9O o /+45°/O°]+0] Laminate 

Wash-out 


wing 

|90°/+45°/0°]-#] Laminate 
Wash -in 


Figure 3 
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OBLIQUE WING ROLL TRIM ANALYSIS 


The main intent of this analysis is to study the performance 
of an oblique wing in roll trim with asymmetric composite 
tailoring. Both cruise and maneuver conditions are explored, the 
cruise case (lg) having a dynamic pressure of 215 psf and an angle 
of attack of -0.25 deg, and the maneuver case (2.25g) having a 
dynamic pressure of 280 psf and a 3 degree angle of attack. Both 
conditions have a Mach number of 0.75. The oblique wing model has 
an aspect ratio of 10 and a taper ratio of 0.4, and incorporates 
the supercritical airfoil OW 70-10-12. The wing structure 
consists of wing skins made of a typical graphite-epoxy composite. 
The composite lay-up and planform shape were shown in figure 3. 
The performance of the wing is measured by four aerodynamic, 
control and structural parameters. Aerodynamically , the pressure 
(induced) drag is noted to see if aeroelastic tailoring results in 
an increase or decrease in drag for the wing. From the controls 
viewpoint, the ability of the ailerons to generate a rolling 
moment (control effectiveness) and the hinge moments on the 
control surfaces are used to measure performance. The hinge 
moments dictate the actuator system for the wing. A decrease in 
hinge moment could result in a lighter actuator system, which is a 
benefit because of a decrease in weight. Structurally, the stress 
level f in the composite skins is noted, defined as 

where o and T are the stresses in the composite layer, and X, Y 
and S are material constants (ref. 12). Before noting how these 
performance parameters are affected by aeroelastic tailoring, 
conditions for roll equilibrium are first obtained by numerous 
aeroelastic computational runs. Figure 4 shows combinations of 
aileron deflection 6 and laminate orientation angle 0 required 
to trim the oblique wing in roll. 
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AERODYNAMIC PERFORMANCE 


Figure 5 plots the pressure drag coefficient versus 
the laminate orientation angle for the oblique wing in roll trim. 
An aileron deflection angle is associated with each laminate 
orientation angle for cruise and maneuver according to figure 4. 
For both cruise and maneuver, the pressure drag remains relatively 
flat. This occurs because the twist distribution across the wing 
is basically the same for the roll-trimmed oblique wing regardless 
of what 6-6 combination is used to achieve that roll 
equilibrium. The drag at 0=20 deg is about 3 or 4 counts higher 
than at 0=0 deg for the cruise and maneuver conditions (one drag 
count equals a drag coefficient of 0.0001). However, it must be 
remembered that the pressure drag does not include boundary-layer 
effects or drag from flow separation. Referring to Figure 4 
again, a fair amount of aileron deflection is required for small 
laminate orientation angles, especially for the maneuver case. We 
would suspect that higher aileron deflections would result in a 
larger boundary layer and a greater likelihood of flow separation, 
which would result in an increase in drag not accounted for in the 
aerodynamic analysis of FL022. Thus, aeroelastic tailoring could 
potentially result in less drag because of the reduction in 
aileron deflection needed for roll trim. 
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CONTROL PERFORMANCE 


The hinge moment coefficient versus laminate orientation 
angle is plotted in Figure 6 for the oblique wing in roll trim. The 
hinge moment would determine the actuator system needed for the 
control surfaces. Since the acutator would be the same for each 
aileron, consider the higher loaded aileron, which is on the aft 
wing. We see that for both cruise and maneuver the hinge moment 
is reduced as the composite laminate is rotated due to the 
reduction in aileron deflection. Because of the reduced hinge 
moment, a smaller, lighter actuator could be used giving a weight 
savings. Aeroelastic tailoring can thus give a performance 
advantage by not only reduced aileron deflection but also a weight 
reduction by the resulting decrease in hinge moments. Additional 
results not shown here also indicate the ailerons will not suffer 
any significant reduction in their ability to produce a rolling 
moment if the wing is aeroelastically tailored. 


HINGE MOMENT COEFFICIENT VERSUS 0 
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STRUCTURAL PERFORMANCE AND CONCLUSIONS 


The effect on the stress level in the composite skins due to 
changing the laminate orientation angle is shown in Figure 7 for 
the oblique wing in roll trim. The figure depicts the maximum 
stress level occurring in the composite laminate, which generally 
occurs in the composite layer whose fiber direction is directed 
mainly along the chord of the wing. It is seen that the maximum 
stress level increases as the composite laminate is rotated. 
This is viewed as a disadvantage because a higher stress level 
would imply that the skin thickness must be increased to obtain 
the desired strength and factor of safety, resulting in more 
weight. 

Thus, performance trade-offs do exist in the application of 
aeroelastic tailoring to the oblique wing. Tailoring the wing 
results in a decrease in the aileron requirements on the oblique 
wing for roll trim, leading to a reduction in aileron hinge 
moments. This implies a weight reduction since a smaller actuator 
could be used. The decreased aileron deflection could also mean 
that aeroelastic tailoring gives a drag reduction because of the 
smaller boundary layer and less likelihood of flow separation 
associated with less aileron deflection. However, aeroelastic 
tailoring also results in an increase in the stress level in the 
composite wing skins, which could result in a weight increase to 
maintain the desired strength. Overall it appears that a 
performance increase is obtained by aeroelastic tailoring. Since 
trade-offs exist, the use of an integrated design approach 
incorporating aerodynamic, structural and control considerations 
would be beneficial (or necessary) for designs with aeroelastic 
tailoring . 
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BACKGROUND 


The elastic semispan wing used in the present study is from the NASA program Drones for 
Aerodynamic and Structural Testing (DAST) and is the right wing panel from the second 
Aeroelastic Research Wing (ARW-2). The DAST ARW-2 wing was designed to flutter within the 
flight envelope of the drone aircraft to which it was to be attached. The wing had an aspect ratio 
of 10.3, a leading-edge sweep angle of 28.8°, and a supercritical airfoil. The wing was designed 
for a cruise condition of M = 0.80, Cl = 0.53 (a = 1.3°), and an altitude of 46,800 feet (q = 

1 27 psf). The wing had three hydraulically driven trailing-edge control surfaces and was 
instrumented with dynamic pressure transducers and accelerometers. 

An unusual transonic instability was encountered near M = 0.9 during an unsteady pressure test 
of the DAST ARW-2 wing in October 1983 in the Langley Transonic Dynamics Tunnel (TDT). 
This unusual transonic instability boundary was predicted using a subcritical response 
technique. This instability was predicted to occur at an almost constant Mach number of 0.9 for 
all dynamic pressures tested. The wing motion was primarily first wing bending mode response 
and was angle of attack dependent. Single degree-of-freedom bending mode oscillations have also 
been encountered during experiments with several other aircraft configurations. These 
oscillations have been observed on a low aspect ratio wing with subsonic airfoil shape, on the B- 
1 A during a wind-up turn, on the canard of the HIMAT aeroelastic model at negative angle of 
attack, and on a forward swept wing force model panel at a negative angle of attack. Linear 
theory flutter analysis was unable to predict any of these oscillations. (Fig. 1.) 


• DAST ARW-2 Wing 

• Aeroelastic wing designed to flutter in flight envelope 

• High-aspect-ratio supercritical wing 

• Design cruise Mach number = 0.8 

• Wing heavily instrumented for unsteady pressure data acquisition 

• Wing tested in TDT, October 1983 (AIAA 85-0598-CP) 

• Instability unexpectedly encountered 

• Boundary predicted using subcritical response technique 

• Motion predominantly 1st wing bending 

• Boundary at almost constant Mach number of 0.90 

• Angle-of-attack sensitive 

• Similar phenomena observed for other aircraft 


Figure 1 
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DAST ARW-2 WING IN TDT 


Figure 2 shows the wing and fuselage configuration mounted in the wind tunnel. The elastic 
semispan wing used in the present study is the DAST ARW-2 right wing panel. A half-body 
fuselage was used to simulate the drone fuselage. This fuselage had shorter nose and tail sections 
than does the drone fuselage since no supersonic tests were to be made. The center section of the 
fuselage was similar to the actual drone fuselage in both diameter and wing location to generate 
the proper airflow over the inboard section of the wing. Both the fuselage and the wing were 
mounted on a remotely controlled turntable mechanism located on the tunnel sidewall. 



Figure 2 
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DAST ARW-2 WING 


The wing planform is shown in figure 3. The wing had an aspect ratio of 10.3 with a leading- 
edge sweep angle of 28.8°. The wing was equipped with three hydraulically driven control 
surfaces, two inboard and one outboard. For this test, the inboard surfaces were held fixed at 0° 
deflection and only the outboard surface was deflected statically. The outboard surface hinge line 
was located at 77 percent of local chord. 

Also shown in figure 3 are the locations of the wing instrumentation. The instrumentation 
consisted of 191 dynamic pressure transducers and 10 accelerometers. In addition, strain 
gauge bridges were located near the wing root to measure bending moments. Small 
potentiometers were used to measure the control surface angular displacement. The model angle 
of attack was measured by a servo accelerometer that was mounted near the wing root. Both 
steady and unsteady pressures were obtained using differential pressure transducers referenced 
to the tunnel's static pressure. Streamwise rows of upper and lower surface pressure orifices 
were located at six span stations: it = 0.274, 0.476, 0.599, 0.707, 0.871 and 0.972. The 
fifth row at r\ = 0.871 lies along the mid-span of the outboard control surface. All of these 
surface orifices were connected to pressure transducers by matched tubes having an inner 
diameter of 0.040 inch and a length of 18 inches. In order to determine the tube transfer 
functions needed to correct the unsteady pressure data from these matched-tube transducers, 
simultaneous measurements were also obtained from a row of in situ transducers mounted on 
the wing upper surface at - 0.875, parallel to the fifth row of surface orifices. Dynamic 
wing deflections were determined using the 10 accelerometers. 



Figure 3 
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INSTABILITY BOUNDARY PREDICTED DURING FIRST TEST 


During the first test of the ARW-2 wing in the TDT an unusual wing instability, with motion 
similar to the wing first bending mode, was encountered. The boundary was determined for a 
wing angle of attack and control surface deflection of 0° as shown in figure 4. Also shown in 
figure 4 as a solid line is the predicted linear theory (doublet lattice theory) flutter boundary, 
which is of a conventional nature. The measured boundary was determined using a familiar 
subcritical response technique known as peak-hold. The boundary was predicted to occur at a 
nearly constant Mach number of 0.90 beginning at a low dynamic pressure of about 50 pounds 
per square foot (psf) (R = 874,000) and rising nearly vertical to over 300 psf (R = 
5,300,000). The observed wing motion during the instability was similar to the wing first 
bending mode, the frequency of which was 8.3 Hz wind-off. The instability frequency was 8.6 
Hz at the lowest dynamic pressure and increased to about 13 Hz at the highest dynamic pressure. 

Because of recent interest in angle-of-attack effects and shock induced effects on wing 
instabilities, several additional test runs were made. These runs include variation of the wing 
angle of attack, comparison using air or Freon as the test medium and comparison with and 
without a transition strip near the wing leading edge. The instability was found to be sensitive 
to variation in angle of attack and, generally, the minimum damping occurred at or near zero 
wing root angle of attack. In figure 4 the solid symbol indicates the Mach number and dynamic 
pressure where the comparison tests were made. The results showed no significant difference 
in the instability boundary for tests in air or Freon. Reynolds number values in Freon are 
approximately 3.1 times greater than those obtained in air. There were also no significant 
differences for tests in Freon with or without a transition strip. 


o =0°, 6 m = 0° 



Figure 4 
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DAST ARW-2 WING RETEST 


A second wing tunnel test was performed on the DAST ARW-2 wing to further investigate the 
unusual instability. The primary purpose of the test was to establish firmly the existence of 
the instability boundary and to gather wing response data and dynamic pressure measurements 
to help understand the mechanism forcing the wing oscillations. A secondary purpose was to 
design an active control system to suppress the wing response using the outboard control 
surface. 

Dynamic pressures and wing deflections were measured for a large number of test conditions in 
the TDT using Freon as a test medium. Data were taken at Mach numbers from 0.5 to 0.96 for 
two stagnation pressures. For a Mach number range of 0.8 - 0.96 the two stagnation pressures 
gave a dynamic pressure variation of 125 - 166 pounds per square foot (psf) and 260 - 340 
psf. These two stagnation pressures will be referred to as the low and high density conditions. 

Wool tufts were used to visualize the flow patterns on the wing in the instability region. Static 
wing tip deflection was measured over the range of test conditions. The effect on dynamic wing 
response of wing angle of attack, static outboard control surface deflection and a spanwise fence 
on the lower surface were investigated. All dynamic wing response data presented is obtained 
from the rear wing tip accelerometer. (Fig. 5.) 


• Purpose 

• Obtain further information on mechanism of non-classical instability 

• Acquire experimental database 

• Assess feasibility of active suppression of wing response (Al AA 87-0881-CP) 

• Method 

• Probe instability region 

• Record unsteady pressure and wing acceleration 

• Flow visualization 

• Measure wing deflection 


Figure 5 
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WING-TIP ACCELEROMETER PEAK-HOLD RESPONSE - 

MACH EFFECTS 


Figure 6 shows the peak-hold results from the wing tip accelerometer for both the low and high 
density conditions. The wing angle of attack and control surface deflection were held at 0°. The 
data show that no instability was found but instead a region of high dynamic wing response was 
observed. For the lower density condition (q = 1 25 - 1 66 psf) the wing motion reaches a 
maximum at M = 0.93 and then rapidly decreases with increasing Mach number. The same 
trend occurs for the higher density condition (q = 260 - 340 psf) with maximum wing motion 
occurring near M = 0.92. The observed wing tip maximum dynamic amplitudes are noted in 
figure 6. At the lower density condition, the amplitude of the wing tip motion was 
approximately 2 inches peak-to-peak. At the higher density condition, which has double the 
dynamic pressure, the amplitude of the wing tip motion doubled to approximately 4 inches 
peak-to-peak. 

Also shown in figure 6 at the higher density condition and M = 0.92 is a single point for a = -1° 
where 6 inches peak-to-peak amplitude of wing tip motion was observed. At this condition the 
wing motion was so severe that the tunnel bypass valves were opened to rapidly reduce the 
dynamic pressure and associated wing motion. 



1 

Peak-hold 

response 



q = 260 - 340 psf 



Figure 6 
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WING-TIP ACCELEROMETER PEAK-HOLD RESPONSE - 
ANGLE OF ATTACK EFFECTS 


Figure 7 shows the peak-hold results from the wing tip accelerometer for three wing angles of 
attack at the lower density condition. The mean control surface deflection was held at 0°. As 
shown in the figure, the maximum dynamic wing response occurred for a wing angle of attack of 
0° at M = 0.93. Changing the wing angle of attack to 2 and -2 degrees decreased the maximum 
wing response and shifted the corresponding Mach number to 0.94. Similar angle of attack 
trend results were seen at the higher density condition. Data were taken for wing angles of 
attack of 1, 0 and -1 degrees up to M = 0.9 and fell within the scatter of the experimental data 
observed at 0° as shown in figure 6. Therefore it is believed that for this configuration the 
maximum wing response occurs when the wing angle of attack is nominally at 0°. 


6 m = 0°, q = 125 - 166 psf 



Figure 7 
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WING-TIP ACCELEROMETER PEAK-HOLD RESPONSE - 
OUTBOARD CONTROL EFFECTS 


Figure 8 shows the peak-hold results from the wing tip accelerometer for three outboard 
control surface mean deflection angles at the lower density condition. The wing angle of attack 
was set at 0°. The figure shows a small increase in wing response for the control surface 
deflection of 6° (trailing-edge down). However, a significant reduction in wing response is 
shown for a deflection of -6°. The wing tip maximum response peak is reduced by half and 
shifted to a lower Mach number of 0.91 . 


a = 0°, q = 125 - 166 psf 



Mach number 


Figure 8 
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DAST ARW-2 WING WITH LOWER SURFACE SPANWISE FENCE 


In an attempt to disturb the flow and change the dynamic wing response, a 1 /2-inch high 
spanwise fence was attached to the lower surface at approximately the 60% local chordline as 
shown in figure 9. The fence ran from the wing planform break (ri = 0.426) to within 5 inches 
of the wing tip (t| = 0.956). The fence was made up of 5 separate one foot-long pieces of 
aluminum placed end to end to minimize increasing the wing stiffness. 



Figure 9 
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WING-TIP ACCELEROMETER PEAK-HOLD RESPONSE - 
SPANWISE FENCE EFFECTS 


The effect of the fence on the wing tip accelerometer peak-hold response at the lower density 
condition is shown in figure 10. The fence has a significant effect upon the wing response, 
lowering the amplitude of maximum wing motion and shifting the peak value to a lower Mach 
number of 0.90. 


a = 0°, 6 m = 0°, q = 125 - 166 psf 



Figure 10 
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DAST ARW-2 WING WITH WOOL TUFTS 


Wool tufts were placed on the upper and lower wing surfaces for several test runs to visualize 
the flow patterns on the wing. The tufts were placed on eight span stations located at ti = .517, 
.558, .635, .671, .761, .816, .905 and .938, as shown in figure 11. The tufts were one inch 
long and on the six inboard span stations were located at ever 1 0% of local chord. On the two 
outboard span stations the tufts were located between 1 0 and 90% chord at every 20% of local 
chord. 



Figure 11 
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SEPARATED FLOW REGIONS INDICATED BY WOOL TUFTS 


Figure 1 2 lists the regions of separated flow on the wing as indicated by the tuft data for Mach 
numbers from 0.85 to 0.96 at the lower density condition. Upper surface flow separation is 
first indicated at M = 0.88. The region of separated flow expands upstream and outboard as 
Mach number increases to 0.94 and then remains constant to M = 0.96. Flow separation on the 
lower surface is initially indicated at M = 0.90. The region of separated flow expands upstream 
and outboard as Mach number increases of 0.94. At M = 0.96 the region of separated flow on the 
lower surface decreases, moving downstream and inboard. 
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MEAN PRESSURE DISTRIBUTIONS 


Figure 13 shows the mean chordwise pressure distribution at the 87.1% span station for nine 
Mach numbers at the lower density condition. The wing angle of attack and outboard mean 
control surface deflection were 0°. As Mach number increases, a shock develops on the upper 
surface at M = 0.85 and becomes quite strong at M = 0.89. The criteria used to determine 
trailing-edge flow separation from mean pressure measurements is the attainment of negative 
pressure coefficients at the 95% chord location. When negative pressures are sustained aft of 
this location, the flow is considered to be separated. Based upon the mean pressure 
distributions shown in figure 13, it appears that flow separation on the upper surface is 
evident at M = 0.92 and is established strongly at M = 0.94. The lower surface develops a 
strong shock at M = 0.92 and the pressure distributions indicate flow separation at M = 0.96. 

Comparing these data to the separated flow regions indicated by wool tufts shown in figure 12 
leads to two conclusions. The first is that the mean pressure data give an incomplete picture of 
the flow separation. The mean pressure data, taken at i\ = 0.871, does not indicate flow 
separation on the upper surface until M = 0.92 while the wool tufts indicate separation in the 
region of the pressure transducers near M = 0.89. Flow separation on the lower surface is not 
indicated by the mean pressures until M = 0.96 while the tufts indicate separation in the region 
at M = 0.92. The second conclusion is that flow separation, as shown by the tuft data, coincides 
with the occurrence of strong shocks on a surface, as shown by the mean pressure data in figure 
13. This flow separation occurs near M = 0.89 on the upper surface and M = 0.92 on the lower 
surface at the 87.1% span station. 


h = 0.871, a = 0°, 6 m = 0°, q = 125 - 166 psf 



Figure 13 
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STATIC WING TIP DEFLECTIONS 


During the test, measurements of the mean wing tip deflection and twist were made using an 
optical cathetometer instrument focused on a straight line drawn on the outboard tip of the wing. 
The results of the wing tip measurements for the lower density condition at a wing angle of 
attack and mean control surface deflection of 0° are shown in figure 14. The wing tip deflection 
and twist increase as Mach number increases up to a maximum near M = 0.85. At higher Mach 
numbers the wing tip deflection and twist values decrease rapidly as the Mach number 
increases. This agrees with the tuft data which show flow separation beginning on the upper 
surface at M = 0.88, causing loss of lift (see figure 13) and the resulting decrease in wing 
deflection and twist. 


a = 0°, 6 m = 0°, q = 100 - 166 psf 


Positive deflection 




Figure 14 
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INSTANTANEOUS PRESSURE DISTRIBUTIONS 


Figure 15 shows the instantaneous chordwise pressure distribution at the 87.1% span station 
for M = 0.92, a = -1° and 8m = 0°. This is the condition at which 6 inch peak-to-peak wing tip 
motion occurred (figure 6). The instantaneous pressure distributions are shown for the 
maximum and minimum vertical wing tip deflection. Based upon the pressure at 95% chord, at 
the maximum wing tip deflection the flow aft of the shock is separated on the upper and lower 
surfaces. The flow is attached on both surfaces when the vertical tip deflection is a minimum. 

This figure points out an important feature of this dynamic motion. At conditions where large 
amplitude dynamic motion is encountered, the trailing-edge flow begins a pattern of separating 
and reattaching on the wing, which coincides with the shock wave motion. As the Mach number 
is increased above 0.92 the flow behind the shock remains separated (see figure 13) and the 
amplitude of the motion rapidly decreases (see figure 6). Thus it appears that the dynamic wing 
response is related to chordwise shock motion in conjunction with shock induced flow separation 
and reattachment on both the upper and lower surfaces. This conclusion is supported further by 
the results obtained when the spanwise fence was attached to the wing lower surface. The fence 
prevented reattachment of the flow on that surface and the maximum wing motion was found to 
be dramatically reduced as shown in figure 1 0. 


H = 0.871, M = 0.92, a = -l°, 6 m =0° 

o □ Maximum deflection j— ? 

• ■ Minimum deflection H 
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PRESSURE VARIATION WITH TIME 


Figure 16 shows a time history of upper and lower surface pressures at the same span station 
and flow conditions as given in figure 15. All pressures are arbitrarily plotted so that they fit 
near each other. However, the last chordwise pressures on both surfaces are plotted with a zero 
reference line. For the last chordwise pressures, the figure clearly shows the separation and 
reattachment of flow in the trailing-edge region of both surfaces as the pressure values 
fluctuate above and below zero. The shock motion, as indicated by large pressure variations, is 
also shown in the figure. For example, the upper surface'shock can be seen to move from in 
front of 68.0% of local chord to behind 74.2%. The lower surface shock moves from in front of 
46.0% to behind 51.3%. At the bottom of the figure the measured wing root bending moment 
time history is plotted for reference. For the observed motion, the wing root bending moment is 
proportional to wing tip displacement, bending maximum for maximum positive (up) wing tip 
position. 

The alternating separation and reattachment of the flow on the upper and lower surfaces 
explains the discrepancy between the mean pressure and wool tuft data. The mean pressure data 
give an average of the pressure values in the trailing-edge region. If, on the average, the flow 
is attached most of the time, the mean pressure distributions will indicate that the flow is 
attached. The mean data give an accurate indication of separation only when the flow remains 
separated most of the time. Another point to note is that while the wool tufts indicate flow 
separation, they are inadequate for indicating the subsequent flow reattachment. 


n =0.871, M = 0. 92, a = -1°, 6 m =0° 
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DEMONSTRATION OF ACTIVE SUPPRESSION OF WING RESPONSE 


The design of a controller to actively suppress the wing response presented several challenges. 
Flutter analysis based upon linear potential flow aerodynamic theories did not predict the wing 
response. The approach adopted was to develop approximate linear models by utilizing forced 
response data taken during the previous tunnel entry in October 1983. Key transfer functions 
were estimated from the forced response data for a range of test points as the region of dynamic 
wing response were approached. Because of the uncertain nature of the transfer function 
estimates upon which the control law was based, additional transfer function and controller 
performance data were gathered during the test and used to modify the controller to obtain 
satisfactory performance. Computer algorithms were written to give near real-time 
assessment of controller performance. 

During the test several sets of fast sine sweep data were taken and averaged for improved 
transfer function estimation. The loop transfer function was estimated using Fast Fourier 
Transform (FFT) techniques and, with the feedback loop open, provided a near real-time 
assessment of controller performance. The controller was modified as necessary and the control 
loop closed. Fast sine sweep data were again taken to evaluate controller performance. (Fig. 
17.) 


• Challenges 

• Mathematical models more uncertain than normal 

• Requirement to estimate key transfer functions from experimental data 

• Development/assembly of algorithms that allow near real-time assessment 
of controller performance during test 

• Method 

• Upgrade transfer function estimates 

• Evaluate controller performance with loop open 

• Confirm modified controller performance by closing loop 


Figure 17 
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CLOSED LOOP BLOCK DIAGRAM 


The inability of flutter analyses based upon linear potential flow aerodynamic theories to 
predict the apparent instability led to a decision to attempt to design a control law based upon 
forced response data from the previous TDT entry (October 1983). Figure 18 depicts the block 
diagram for the closed-loop single-input/single-output system. The transfer function 
estimates, z/8a, were obtained from response data due to fast sine sweep inputs into the aileron 
actuator, Ga. These data were available for a number of test points in the range of Mach number 
0.7 < M < 0.85. Unfortunately, the accuracy of the estimates became more suspect as the 
apparent stability boundary, where good estimates were most needed, was approached. 
Nevertheless, these transfer function estimates were employed in the preliminary design of a 
control law, H, to add damping to the critical mode. 



Figure 18 
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ROBUSTNESS MAXIMIZATION 


The objective of the control law design activity was to add damping in the critical frequency 
range without degrading the stability characteristics of the higher frequency modes. The 
approach taken was to define, for a range of test points, the compensation, H(s), required to 
modify the amplitude and phase of the accelerometer outputs such that damping was added in a 
frequency band centered at the frequency of the first elastic mode. The control law form chosen 
is shown in figure 19. 

T(s) is a fixed low pass filter chosen to attenuate the feedback control so that higher frequency 
elastic modes were not affected. The factor adjacent to T(s) concentrates the control effort in 
the frequency region, an, of the first elastic mode; the frequency, coi(M,q), depended upon 
Mach number, M, and dynamic pressure q. The parameter £ was fixed at 0.2. The choice, 
£=0.2, was made to confine the control activity in a narrow band around coi . In retrospect, a 
larger value of £ would probably have resulted in more attenuation of the response due to 
turbulence and could have been used safely since the frequency of the second elastic mode was 
about 32 Hz. The remaining part of H(s) factor has variable coefficients (Di(M,q), i=1 ,5} 
which were used as control design variables to allow proper amplitude and phasing for 
robustness. 

Values for the variable parameters were found for each test point such that the minimum 
singular value of the return difference transfer function was maximized subject to gain and 
phase margin constraints. The search for the optimizing set {Di* (M,q)} was done using a 
nongradient constrained optimization approach. 


Find values for the design variables D . where 

(S 2 + D, S + DJ 

H=D 1 

(S + D 4 S + D $ ) 

V v / 

Phasing adjustment 
Such that 

• Minimum singular value is maximized 

• Good gain and phase margins 


tfjS 

S + 2 £ j 

' ^ ' 



T (S) 



Figure 19 
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ACCELEROMETER PEAK-HOLD RESPONSE 


Closed-loop performance is shown in figure 20 in terms of frequency domain peak-hold 
responses. Peak-hold is a subcritical response technique used in predicting flutter boundaries 
wherein the autospectrum of an output is obtained for a block of data. Subsequent autospectra 
are taken for a number of blocks of data and the peak value at each frequency out of the entire set 
of blocks is retained (held). As the point of neutral stability of a mode is approached, the 
amplitude at its resonant frequency approaches infinity. Thus, by observing the variation in 
the peaks held, or typically, their reciprocals, as a function of test condition changes, one 
obtains an indication of changes in damping ratio and a prediction of where an instability might 
occur. 

The figure shows peak-hold responses with the control system loop open and closed both for 
turbulence only and turbulence plus fast sine sweep excitation. It is seen that the controlled 
cases have lower responses than the uncontrolled cases over the entire Mach number range. The 
controlled response is sharply attenuated at the first elastic mode resonant frequency and 
somewhat amplified on either side of it as compared with the uncontrolled case. The 
amplification at the slightly higher frequency is due to a lightly damped mode introduced by the 
controller. Increasing the controller badwidth by selection of a lager value for £, as mentioned 
earlier, would probably have allowed further reduction in the controlled response. Control 
effort due to turbulence excitation varied from a peak feedback signal magnitude of 
approximately 0.25° at M = 0.70 to 2.3° at M * 0.92 for these cases. 



q = 100 to 151 psf 


Figure 20 
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CONCLUDING REMARKS 


The DAST ARW-2 wing had been tested previously in the NASA Langley TDT and an unusual 
instability boundary was predicted based upon subcritical response data. Contrary to the 
predictions, no instability was found during the present test. Instead a region of high dynamic 
wing response was observed which reached a maximum value between Mach numbers 0.92 and 
0.93. The amplitude of the dynamic response increased directly with dynamic pressure. 

The response appears to be related to chordwise shock movement in conjunction with flow 
separation and reattachment on the upper and lower wing surfaces. The onset of flow separation 
coincided with the occurrence of strong shocks on a surface. Instantaneous pressure 
distributions indicated that the flow was intermittently separating and reattaching near the 
trailing edge under conditions of maximum wing motion. The dynamic wing response was 
sensitive to angle of attack, with maximum motion occurring near a = 0°. Static deflection of 
the outboard control surface significantly decreased the dynamic response for 8m = -6°. A 
spanwise fence installed on the lower surface to disturb the flow pattern resulted in a 
significant decrease in dynamic wing response. 

A controller was designed to suppress the wing response. The control law attenuated the 
response as compared with the uncontrolled case and added a small but significant amount of 
damping from M = 0.70 to M = 0.92 for the lower density condition (q = 100 - 151 psf). 

The unsteady pressure and response data acquired during this test constitute a valuable data base 
to be used for further study of this unusual phenomena and for validation of unsteady CFD codes. 
(Fig. 21.) 


• No "hard flutter" point obtained although amplitude increases as 
dynamic pressure increases 

• Response coincides with onset of shock-induced separation and 
reattachment 

• Response affected by angle-of-attack, outboard control surface 
position, and a spanwise fence on the lower surface 

• Active suppression successfully demonstrated, significant damping 
added with feedback loop closed 

• Unsteady pressure data acquired during dynamic response for 
further analysis and use in code validation work 


Figure 21 
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THE ROLE OF SHOCK INDUCED TRAILING EDGE 
SEPARATION IN LIMIT CYCLE OSCILLATIONS 


In figure I Limit Cycle Oscillations (LCO) are defined as limited 
amplitude oscillations which are self-sustaining and are produced by a 
structural/aerodynamic interaction. The role of shock induced trailing 
edge separation (SITES) in this phenomenon, as will be discussed in 
this paper, is to act as a non-linear spring which triggers and drives 
the LCO. The appearance of SITES coincides with the classical trailing 
edge pressure divergence which is a well-known indicator of buffet 
onset. Because of this, LCO has previously been referred to as buffet 
or, at transonic speeds, Mach buffet. The conditions for transonic LCO 
are moderate incidence, usually less than 10°, and high transonic Mach 
numbers, ranging from 0.8 to 1.1. Lowly damped vibration modes tend to 
respond provided they have the proper characteristics to couple with 
the SITES type flow. These conditions frequently occur near flutter 
boundaries which creates considerable anxiety in both test pilots and 
engineers, especially when the response is in a single mode with a 
distinct frequency. Since amplitudes can become quite large, even 
though they are limited, the knowledge of LCO boundaries becomes very 
important for efficient flight flutter testing. 


• LCO Is Defined As A Limited Amplitude Oscillation Which Is Self- 
Sustaining And Results From A Structural /Aerodynamic Interaction 

• Shock Induced Trailing Edge Separation (SITES) Coincides With 
The Classical Trailing Edge Pressure Divergence 

• The Conditions Are Moderate Incidence And High Transonic 
Mach Numbers 

• Lowly Damped Vibration Modes Tend To Respond If They Have 
The Proper Characteristics 

• Response Is Single Mode With Distinct Frequency And Limited 
Amplitude 


Figure 1 
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WING BENDING AND TORSION MODE RESPONSE FROM 
THE F-lll TACT BUFFET ANALYSIS AT M=0.8, A LE =26 DEG 


The results shown in figure 2 represent the comparison of 
predictions and flight test results for buffet response of the F-lll 
tact aircraft at M=0.8 and a leading edge sweep of 26° (ref. 1). The 
results are wing tip RMS accelerations for the first symmetric wing 
bending mode and the first wing torsion mode group which includes both 
asymmetric as well as symmetric modes. Several altitudes are 
represented which range from 13K ft to 28K ft, for angles varying from 
7° to 12°. The interesting feature to note in the comparisons is that 
(1) the bending mode responses are very well predicted in terms of both 
angle and altitude effects whereas (2) the torsion mode responses are 
consistently underpredicted with exception of the lowest altitude. 

Since the prediction method did not allow for a coupling between the 
buffeting flow and the structural response, these results led to an 
investigation of what kind of coupling could exist with the torsion 
modes. Thus, the role of SITES in LCO was developed as an answer to 
this question. 


(Reference 1) 



Angle Of Attack, Deg. Angle Of Attack, Deg. 

First Wing Bending Wing Torsion Mode 

Mode Response Response 

Figure 2 
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STATIC PRESSURE CHANGES ON THE MODEL AND AIRPLANE 
WING UPPER SURFACE FOR TRANSITION TO SITES 


A clue to a possible means for the coupling can be deduced from 
the static pressures in figure 3. The occurrence of trailing edge 
pressure divergence at about 10° also corresponds to a large forward 
movement of the upper surface main shock as shown by comparing the 
pressures in figure 3 at 9° and at 10°. It will also be noted that the 
forward shock movement for the aircraft is much larger than that for 
the model. The condition of this transition is the occurrence of 
shock-induced trailing edge separation (SITES) which was extensively 
discussed by Cunningham, et al in reference 2. It was shown that this 
transition was accompanied by a step change in pitching moment with 
either increasing or decreasing angle-of-attack. With increasing 
angle, the forward shock movement produced a loss in lift forward and 
the trailing edge divergence produced a gain in lift aft. The net 
result was to provide a step change in pitching moment that was nose 
down for increasing angle. For decreasing angle, the opposite took 
place and produced a step change in pitching moment that was nose up. 


(Reference 1) 



1/6 Scale 

Model 

r MAC = 14 xIO 6 

FLT Test 

r MAC = 25 x 10® 

M = 0.8 
Ale = 26 De 9- 



a = 9.1 Deg 
Pre-SITES 


CX = 10 Deg 
Post-SITES 


Figure 3 
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EFFECT OF MEAN ANGLE, FREQUENCY AND AMPLITUDE 
ON MEASURED UNSTEADY PRESSURES 


Dynamic unsteady investigations of SITES flows were conducted by 
Triebstein (ref. 3) for an NACA 0012 airfoil oscillating in pitch at 
various mean angles, frequencies and amplitudes as shown in figure 4. 
These results are the upper surface unsteady pressure distributions for 
the first harmonic. The effect of varying mean angle is quite 
pronounced where at a = 0° the distributions are of conventional 
transonic characteristics, but at a = 5° the distributions are totally 
different. The most notable change is the trend toward a more nose- 
down pitching moment. The influence of increasing frequency simply 
shows an expected increasing lag at either angle. Increasing amplitude 
at a = 5° with SITES shows a decrease of the shock motion peak amplitude 
when normalized by the pitching amplitude. This demonstrates the 
tendency to maintain a constant incremental force with increasing pitch 
increments so that a limited amplitude motion is inevitable. That is, 
the viscous damping which is a function of amplitude increases until it 
balances the destabilizing fixed incremental force due to the SITES 
transition. 


(Reference 2) 

Upper Surface Results For An Oscillating NACA 0012 Airfoil, M = 0.78 



□ a = 0 ° 

k = 0.073 
R e = 2.6x10® 

o or = 0° 
k = 0.024 
R e = 2.6X10®J 


M = 0.78 A« = 1.0° 


a = 5° 
k = 0.024 
R e = 1.1x10® 

a - 5° 
k = 0.073 
R e = 1.1 x10®J 


WITH 

SITES 


Varying Mean Angle 
And Frequency 
With Fixed Amplitude 


-10 


NO 

SITES 



M = 0.78 
k = 0.072 
a = 5° 

R e = 1.1X10® 
o A a = 0.5° 

□ A«= 10° 

■ A cr= 1.5° 

• A « = 2.5° 


Varying Amplitude 
At Mean Angle Of 5 Deg. 
With SITES 


Figure 4 
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SCHEMATIC DIAGRAM OF THE VARIATION OF THE INCREMENTAL PITCHING 
MOMENT DUE TO SITES WITH ANGLE AND PITCH RATE 


The step change in pitching moment just described in figures 3 and 
4 can be cast in the form of a non linear spring that provides a step 
increase in resisting stiffness as the airfoil passes through SITES. 
This form is illustrated in figure 5 for a hypothetical situation in 
which, for simplicity, only the step change is present. (Normally, 
the pitching moment variation with angle is a nearly constant slope 
with the discontinuity superimposed to produce a shift in 
characteristic at the SITES transition point.) The influence of pitch 
rate that provides the hysteresis is also illustrated. For positive 
pitch rate, SITES is delayed to a higher angle, whereas for negative 
pitch rate, re-attachment is delayed to a lower angle. 


A PM (a) 
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HYPOTHETICAL TIME HISTORY OF AN AIRFOIL 
UNDERGOING SITES INDUCED TORSIONAL LCO 


The non-linear spring described above can potentially produce a 
limit amplitude, self-sustaining oscillation. How this is possible can 
be described by considering an airfoil with a torsion spring undergoing 
a maneuver of increasing angle until it encounters SITES as shown in 
figure 6. Slowly increasing incidence at angles below that of SITES 
allows the torsion spring to attain a continuous state of equilibrium 
with aerodynamic pitching moment. When SITES is reached, a sudden 
nose-down increment is imposed on the aerodynamic pitching moment which 
will tend to reduce wing incidence. This will be a dynamic negative 
pitch rate which will delay re-attachment and permit the nose-down 
moment to put work into the system. At some point, re-attachment does 
take place and the nose-down moment disappears. Accelerations become 
negative and the wing experiences a reduced pitch rate until it reaches 
zero and begins nose-up motion. Positive pitch rate now takes over 
which will produce a delay in SITES and allows an overshoot of the 
initial starting point due to stored elastic energy during the down 
stroke. When SITES does occur, the cycle then repeats itself. 
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A MATH MODEL FOR LCO 


The math model for calculating LCO response was developed and 
presented in reference 1. The governing equation for a non-linear 
single DOF spring-mass system is shown in figure 7. The LHS of the 
equation is the conventional linear equation form for a generalized 
coordinate response. The RHS is composed of the non-linear step force 
function as defined for the vibration mode of interest. This force is 
obtained by integrating the incremental pressure changes due to SITES 
transition with the mode shape to produce a generalized force. The 
boundary conditions require that at time zero, the system is in static 
elastic equilibrium (zero velocity) and that the normalized 
displacement is equal to zero at the static transition point for SITES. 
That is, the step function change on the RHS is suddenly imposed on the 
system at time zero. 


• The Non-Linear Single DOF Spring-Mass System 
Mj q + 2 Mj Wj 6j q + Mj u)j 2 q = AFj (q, q) 


• Where 

Mj, Uj, 8j, q = Mass, Frequency, Damping and Response 

For The ith Mode 

A Fj(q, q) = Non-Linear Incremental Force 
(Due To SITES) 


• Is Subject To The Initial Conditions 
q (t = o) = qtrans (q = o) 
q (t = o) = o 

Figure 7 
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SCHEMATIC DIAGRAM OF THE VARIATION 
OF AF^q,^) NEAR SITES 


The plot shown in figure 8 illustrates that the form of variation 
of the step change in generalized force is identical to that shown for 
pitching > moment in figure 5. The influence of pitch rate in the 
hysteresis is also identical for wing modes that have significant 
torsional motion. This is consistent with observations to date which 
indicate that LCO occurrences which are coincident with SITES tend to 
concentrate on torsional modes or modes with some pitching content. 


AFj (q, q) 
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A MATH MODEL FOR LCO (contd.) 


A finite difference solution was developed in reference 1 for the 
equation of motion discussed in figure 7. The form of the solution is 
shown in figure 9 where the parameters have all been cast in non- 
dimensional form with exception of € and q n . These terms have the 
dimension of length. The two equations are applicable depending on the 
response value relative to SITES transition. Hysteresis is included in 
the model through the use of a transition point that is sensitive to 
wing motion velocity, q. The initial conditions are defined such that 
the transition actually occurs at the third time step, q 2 * 


• A Finite Difference Solution To The Non-Linear Equation Is 

9 n+1 = 1 + ^g [ A 2 e + (2 - A 2 )q n “ 0 - AS )q n-l], 9 n > q Trans (d) 

9 n+1 = TTxS ^ — n “ ^ — AS )q n-l]> 9 n < q Trans (d) 

• Where 

A = At wj e = - 

Mi coj2 

5 = 8j 

A t = Time Step Size 

q n = Response Of The ith Mode At The nth Time Step 

• Subject To The Initial Conditions 

q ° =qi = 0 ’ q2= TTTi 

Figure 9 
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SAMPLE OF LCO CALCULATIONS FOR THE F-lll TACT 


The algorithm just described was programmed and input data were 
developed for the right wing torsion mode of the F-lll TACT airplane 
(ref 1). The conditions were M=0.8 and 2 IK ft altitude. The 
generalized mass, frequency and damping data were obtained from a 
conventional dynamic response and flutter analysis. The estimate for 
AF^ was developed by using the pressure distributions at 9° and 10° to 
obtain the incremental loads due to transition to SITES. The load 
distribution was combined with the mode deflections, h^ (x,y) to 
produce the generalized step force value of € . Estimates for Aq, the 
hysteresis parameter (see fig. 8) , were not possible based on available 
data, hence, a parameter study was conducted by letting Aq vary as a 
fraction of € . A sample plot for Aq = € shown in figure 10 
illustrates how the transient solution quickly approaches a limit 
amplitude oscillatory motion. It is also interesting to note that the 
apparent frequency is slightly higher than the right wing torsion mode 
frequency. (Since the number of time steps shown in the plot 
represents ten cycles, the apparent frequency can be determined by 
simply counting cycles.) 


• The Parameters For The Right Wing Torsion Mode Were 


Mi 


CO 


Si 

e 


l = 


= 565.9 Lbs 

= (14.17 Hz) X 2ir 

= 0.07 (From Flutter Solution) 

— - — f hj (x, y) [Cp 10 (x, y) - Cp 9 (x, y)] d x d y 
M i"i 2 'Awing 
-0.0127 Ft 



Aq = 0.01270} 


A 



0 . 200 . 400 . 600 . 800 . 1000 . 1200 . 

Nt — Number of Time Steps 


Figure 10 
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RESULTS OF LCO CALCULATIONS FOR THE F-lll TACT AIRCRAFT 
AND 1/6-SCALE MODEL RIGHT WING TORSION MODE 


The results shown in figure 11 are for the nominal conditions of 
M=0.8, a=9°-10° as were presented in reference 1. The flight test data 
were extracted from figure 2 by subtracting the predicted response from 
the flight measured response for the wing torsion modes. Since the LCO 
attributed response was assumed to be uncorrelated with the broad-band 
buffet response, this subtraction was done on an RMS basis. The 
calculated full scale data were obtained with the LCO algorithm by 
letting Aq vary as a fraction of € from 0.2 up to 2.0. It is 
interesting to note that up to a ratio of 1.0, the LCO response 
increases and the frequency drops. Within this range, the calculated 
amplitude is within reasonable agreement with the flight data 
considering that the method should be conservative as it currently 
exists. For Aq =2£, the response is zero, which states that the 
conditions for sustained oscillation require that the transition must 
occur in both directions during the cycles; otherwise, the response is 
simply a static displacement. Finally, the 1/6-scale model results 
show that the calculations verify the observation that the model did 
not experience LCO in the wind tunnel test. 


Conf. 

Method 

£(Ft) 

|e (Deg.)| | Aq/e| 

APP. f, Hz 

RWAC. 

g’s Response 

F.S.A.C. 

FLT Test 

- 

- 

- 

- 

— 

15.0 

0.8 

F.S.A.C. 

Calc. 

-0.0127 

0.073 

0.2 

18.4 

0.71 









0.4 

17.0 

1.42 









0.6 

16.4 

1.79 









0.8 

15.1 

2.05 



} 


1 


1 


1.0 

2.0 

14.9 

14.2 

2.34 

0 

1/6-S.M. 

Wind Tun. 

- 

- 

- 

- 

— 

156.0 

a 0 

1/6-S.M. 

Calc. 

-0.000253 

0.0087 

8.4“ 

156.0 

0 


*LCO Response From Testing Is Estimated As The Amount That Exceeds Expected 
Linear Buffet Response 

** Aq Assumed As A Constant Angle Determined By F.S.A.C. For | Aq/e| = 1.0 

Figure 11 
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CONCLUDING REMARKS 


As a result of this investigation, the potential role of shock 
induced trailing edge separation in limit cycle oscillations was 
established. It was shown that the flip-flop characteristic of 
transition to and from SITES as well as its hysteresis could couple 
with wing modes with torsional motion and low damping. This connection 
led to the formulation of a very simple non-linear math model using the 
linear equations of motion with a non-linear step forcing function with 
hysteresis. A finite difference solution with time was developed and 
calculations were made for the F-lll TACT airplane. Static pressure 
data for the F-lll TACT were used to determine the step forcing 
function due to SITES transition. However, no data were available for 
the hysteresis hence a parameter study was conducted allowing the 
hysteresis effect to vary. Very small hysteresis effects, which were 
within expected bounds, were required to obtain reasonable response 
levels that essentially agreed with flight test results. Also in 
agreement with wind tunnel tests, LCO calculations for the 1/6-scale 
F-lll TACT model showed that the model should not have experienced LCO. 


• The Role Of Shock Induced Trailing Edge Separation (SITES) In Limit 
Cycle Oscillations (LCO) Has Been Described 

• The Flip-Flop Characteristic Of Transition To SITES and Its Hysteresis 
Has Been Shown To Be A Potential Source Of LCO 

• A Very Simple Non-Linear Math Model Was Assembled and Solved 
With A Finite-Difference Approach 

• The Math Model Used Static Pressure Model Data and Standard Flutter 
Solution Results; However, Hysteresis Data Were Not Available and Had 
To Be Assumed For This Study 

• Calculations For The F-111 Tact Aircraft Showed Order Of Magnitude 
Agreement With Flight Test 

• Calculations Also Showed That The 1/6-Scale Model Should Not Have 
Experienced LCO-Which Agrees With Observations 

Figure 12 
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INTRODUCTION 


The purpose of this presentation is to present a brief status report on the initial application of 
the CAP-TSD computer program for wing flutter analysis. The CAP-TSD program (for 
Computational Aeroelasticity Erogram - Iransonic Small disturbance) is based on an 
approximate factorization (AF) algorithm that is stable and efficient on supercomputers with 
vector arithmetic. The program has been described by Batina* in an earlier presentation. CAP- 
TSD has been used to calculate steady and unsteady pressures on wings and configurations at 
subsonic, transonic, and supersonic Mach numbers. Comparisons of these results with other 
methods and with experimental data have been favorable. However, the CAP-TSD code has been 
developed primarily for aeroelastic analysis. The present paper reports on the initial efforts 
for validation of the aeroelastic analysis capability. The initial applications include two series 
of symmetric, planar wing planforms. Well-defined modal properties are available for these 
wings; this is vital for accurate flutter calculations. In addition, transonic flutter boundaries 
are available for evaluation of the transonic capabilities of CAP-TSD. Additional comparisons 
are also being made with linear theory and with the 2-D code XTRAN2L. (Fig. 1 .) 

*Batina et al., NASA CP- 3022, 1989, Paper No. 4, pp. 63-96. 


• COMPUTATIONAL AEROELASTICITY PROGRAM - TRANSONIC SMALL 
DISTURBANCE 

• PREVIOUS EMPHASIS HAS BEEN ON PRESSURES 

- GENERALLY GOOD RESULTS 

- HAVE CONSIDERED STEADY AND UNSTEADY CASES 

- CONFIGURATIONS 

• PROGRESS REPORT ON AEROELASTIC VALIDATION 

• CONSIDERING SYMMETRIC PLANAR WINGS 

- WELL-DEFINED MODAL PROPERTIES 

- TRANSONIC FLUTTER BOUNDARIES 

• COMPARISONS WITH 2-D CODE XTRAN2L ARE ALSO UNDER WAY (NOT 
PRESENTED) 


Figure 1 
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WINGS 




Two series of wing planforms are being used for the initial flutter calculations with CAP-TSD. 
The first set of wings is a series of swept and tapered wings that are being considered as an 
AGARD standard configuration for aeroelastic analysis. These wings are swept back 45° at the 
quarter chord. They are described further by Dr. E. C. Yates' presentation* of this workshop. 
The wings and test data are presented in NASA TN D-1616, dated March 1963. 

The other wing planform is a clipped delta wing that was used in some early flutter suppression 
studies. It is described in NASA TN D-7544, June 1974, and NASA TR R-450, December 
1975. The leading edge sweep for this wing is 50.5° and it is highly tapered. (Fig. 2.) 

* Yates, E. C., NASA CP- 3.022, 1989, Paper No. 12, pp. 243-26.0. 


• 450 SWEPT WINGS - NASA TN D-1616 AND YATES’ PRESENTATION 


• CLIPPED DELTA WING - NASA TN D-7544 AND TR R-450 


Figure 2 
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The planview of the 45° swept wing is shown in figure 3. The wings were semispan, wind- 
tunnel-wall-mounted models that had a quarter chord sweep of 45° (leading-edge sweep of 
46.3°), a panel aspect ratio of 1.65, and a taper ratio of 0.66. The wings had an NACA 65A004 
airfoil section and were constructed of laminated mahogany. In order to obtain flutter for a wide 
range of Mach number and density conditions, some of the wings had holes drilled through the 
wing to reduce the stiffness. To maintain the airfoil shape, the holes were filled with a rigid 
foam plastic as can be seen in figure 3. 



Figure 3 
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45° SWEPT WING IN THE NASA TRANSONIC DYNAMICS TUNNEL 


One of the 45° wings is shown mounted in the Transonic Dynamics Tunnel (TDT) at NASA 
Langley Research Center in figure 4. The models were tested in air and Freon+ test media. The 
semispan of most of these models was 2.50 feet, which is small compared to the 1 6-foot test 
section of the TDT. The models were tested at zero angle of attack. 


+Freon: Registered trademark of E. I. duPont de Nemours and Co., Inc. 



______ 




Figure 4 
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OBLIQUE PROJECTIONS OF NATURAL VIBRATION MODES 

45° WING WEAK3 


The vibration mode shapes for the 45° wings were not measured, but node lines, frequencies, 
and stiffnesses are available from the report (NASA TN D-1616, 1963). Mode shapes were 
calculated with a finite-element analysis and the wing properties were adjusted to match the 
measured nodes lines and frequencies for the lower vibration modes. Oblique projections for the 
first four modes for wing WEAK3 are shown in figure 5. The modes numbered 1 through 4 
represent first bending, first torsion, second bending, and second torsion, respectively. The 
modal frequencies range from 9.60 Hz for the first bending mode to 91.54 Hz for the second 
torsion mode. Similar mode shapes have been calculated for the wing SOLID2. 



(b) Mode 2, f 2 =38.17 Hz 



(d) Mode 4, f 4 =91.54 Hz 


Figure 5 
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AEROELASTIC TRANSIENT AND LEAST-SQUARES CURVE FIT 


For aeroelastic analysis, the steady-state flow field is first calculated to account for wing 
thickness, camber, and mean angle of attack. The wing is then disturbed with an initial 
condition and free decay transients are calculated. The resulting transients are then analyzed to 
determine growth or decay for aeroelastic stability. Dynamic pressure is changed, and the 
transients computed again to determine the variation of stability with dynamic pressure. 

An example transient for the 45° wing calculated by CAP-TSD is shown in figure 6. All four 
modes used in the analysis were excited by specifying an initial condition for each modal 
velocity which produces a complicated decay record. This record is analyzed using a least- 
squares curve-fit of the response data with complex exponential functions. The program 
utilized is a derivative of the one described by Bennett and Desmarais in NASA SP-415, May 
1975. 


Amplitude 



Figure 6 



COMPONENT MODES FROM CURVE FIT 


The components of the transient presented in the previous figure are shown in figure 7 to the 
same scale as that used in figure 6. The free decay properties of each mode for this condition are 
readily apparent after the least-squares fit, whereas the complexity of the complete decay 
record is such that the stability is not recognizable in the previous figure. The instability of 
the first mode might have been missed unless many more time steps were run. A post- 
processing program of this type is essential to efficient use of these types of programs where 
large resources are used for the CFO flow field calculations. 



Figure 7 
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EXAMPLE OF ROOT LOCUS FROM CAP-TSD RESULTS 
45° WING, M = 0.499 


The potential of this methodology to produce complete root loci for the aeroelastic system is 
illustrated in figure 8. The variation of frequency and damping for all four modes used in the 
analysis is deduced for various dynamic pressures as shown in the figure (note change in 
frequency scale for the higher modes). It is apparent that the first mode increases rapidly in 
frequency and flutters, whereas the damping in the other three modes increases rapidly with 
dynamic pressure. 



Tl, Hz 


Figure 8 



PRELIMINARY FLUTTER CALCULATIONS FOR 45° SWEPT WINGS 


Preliminary flutter calculations for the 45° wing WEAK3 in air are shown in figure 9. The 
circles indicate the measured flutter points which are faired by the solid line. The bottom of the 
dip near Mach 1 .0 is estimated from the no flutter data obtained while going to the point at M = 

1 .07. The squares indicate the results from subsonic kernel function linear theory (program 
FAST). There is very good agreement of the linear theory with the four data points shown, even 
for the point near M = 0.95. The two diamonds indicate two subsonic points calculated using 
CAP-TSD. The two points are in fair agreement with the data. Effort to extend these results to 
other Mach numbers and to obtain direct comparisons with linear theory is continuing. These 
initial results are encouraging however. 


Experiment 



M 


Figure 9 
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CLIPPED DELTA WING IN THE NASA LANGLEY 
TRANSONIC DYNAMICS TUNNEL 


The second wing to be analyzed is a clipped delta wing model that was also tested in Freon in the 
Langley Transonic Dynamics Tunnel. A view of the model mounted in the TDT is shown in figure 
10. The wing has a leading-edge sweep of 50.5°, a panel aspect ratio of 1.24, and a taper ratio 
of 0.142. The airfoil section is a circular arc with a maximum thickness of 0.03. The wing 
was constructed of a load-carrying plate structure with cutouts to simulate a beam structure 
and was covered with balsa wood which was contoured to the required airfoil shape. The model 
also had two slender underwing bodies to simulate engine nacelles. The total mass of these bodies 
was about the same as the total mass of the wing. A fuselage fairing was used to ensure that the 
wing root was outside the tunnel wall boundary layer. Nine natural vibration modes and their 
associated generalized masses were measured for this wing. 



Figure 10 
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CLIPPED DELTA WING FLUTTER BOUNDARY 
NASA TN D-7544 


The experimental flutter boundary for the clipped delta wing is shown in figure 11. 

Calculations with CAP-TSD are under way, but results have not yet been obtained. The figure 
shows a composite boundary obtained by normalizing the boundaries for three wings of differing 
sizes in terms of the flutter speed index. The data for the wing of the previous figure are shown 
as the diamond symbols. The flutter boundary has a significant transonic dip with a minimum 
near M = 0.92, a rapid rise after the dip, and a supersonic level near that of subsonic speeds. 



Figure 11 
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CONCLUDING REMARKS 


In an effort to assess the accuracy of the CAP-TSD program for aeroelastic applications, flutter 
calculations are under way for several wings of two different planforms varying in sweep and 
taper and with thin airfoil sections. One planform is a series of 45° swept wings which have 
been proposed as an AGARD standard configuration for aeroelastic analysis. The other planform 
is a clipped delta wing that was used in some early active controls work. The physical 
properties and experimental flutter boundaries for these wings are well defined for validation 
purposes. Some initial results have been obtained and are encouraging. Further effort to extend 
and refine the results is under way. (Fig. 12.) 


• CAP-TSD AEROELASTIC VALIDATION UNDER WAY 

• HAVE INITIAL RESULTS WHICH ARE BEING EVALUATED AND REFINED 

Figure 12 
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UNSTEADY EULER EQUATIONS 


In transonic flutter problems where shock motion plays an important part, 
it is believed that accurate predictions of the flutter boundaries will require 
the use of codes based on the Euler equations. Only Euler codes can obtain the 
correct shock location and shock strength, and the crucially important shock 
excursion amplitude and^ phase lag. (For a discussion of the importance of 
shocks in transonic flutter, see Ref. 1.) The present study is based on the 
finite volume scheme developed by Jameson and Venkatakrishnan (Refs. 2,3) for 
the two-dimensional unsteady Euler equations. The equations are solved in 
integral form on a moving mesh, Eqs . (1-2). Here the variables p, p, u, v and 
e are the pressure, density, cartesian velocity components, and total energy, 
respectively, and and y^ are the velocity components of the moving boun- 
dary 3ft of an element ft , By applying Eq. (1) to each element or cell 
(i,j), a system of ordinary differential equations is obtained, Eqs. (3), where 
Sjj is the cell area, Q^j is the net flux out of the cell, and Dy repre- 
sents dissipative terms added to damp numerical oscillations (see Refs. 3,4). A 
five-stage Runge-Kutta scheme is used to integrate Eqs. (3) forward in time. 


FINITE VOLUME FORMULATION 


• INTEGRAL FORM ON A MOVING MESH 


f- \\ Wdxdy + J (fdy - gdx) =0 (1) 
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TYPICAL SECTION MODEL 


The wing is modeled as a typical section, with two degrees of freedom 
(bending h and torsion a), as illustrated in Fig. 1. The usefulness of this 
model in capturing the fundamental features of bending-torsion flutter is by now 
well established. In the usual notation, the equations of motion are of the 
form given by Eqs. (5) and (6), where the lift and moment coefficients Cl and 
depend on the motion of the airfoil. Because we will consider finite 
(rather than infinitesimal) amplitude motion, the superposition principle cannot 
be used. In the present study, Cl and C^ are calculated numerically from 
the unsteady pressure coefficient on the airfoil surface at the end of each time 
interval, obtained from the numerical solution of the Euler equations. It 
should be emphasized that the equations of motion are nonlinear through the 
dependence of Cl and C^ on the motion h,a of the airfoil (and its time 
history) . 
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METHOD OF SOLUTION 

Aeroelastic stability is determined by integrating the equations of motion 
for the coupled fluid-structure system. The structural equations are first 
transformed to normal coordinates r) r • Eqs • (7)-(8), where the columns of [$] 
are the eigenvectors of the free vibration problem. The structural integrator 
is based on the convolution integral solution, Eq. (9), and the generalized 
aerodynamic forces Q r are assumed to vary linearly within each time step At . 
Because the multi-stage Runge-Kutta scheme used to integrate the unsteady Euler 
equations was found to be sensitive to the manner in which the airfoil boundary 
condition was updated and the mesh moved, the structural integrator has been 
imbedded within the Runge-Kutta scheme in the Euler code. This permits an ef- 
ficient implementation of the exact airfoil boundary condition, Eq. (10), on 
the instantaneous position of the airfoil, given by B(x,y,t) = 0. Nonref lective 
boundary conditions are used in the far field. 


• COUPLED EQUATIONS FOR FLUID & STRUCTURE ARE 

INTEGRATED NUMERICALLY USING NORMAL COORDINATES 


{q} = [4>]{n} 


(7) 
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STRUCTURAL INTEGRATOR IS BASED ON CONVOLUTION INTEGRAL SOLUTION 
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• STRUCTURAL INTEGRATOR IS IMBEDDED IN FIVE-STAGE 
RUNGE-KUTTA SCHEME FOR EULER EQUATIONS 

- EXACT AIRFOIL B.C. IS SATISFIED 

£7 = 0 or I 7 + u»VB = 0 

Dt ot 

on B(x,y,t) = 0 

- MESH IS MOVED AT EACH STEP 
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NUMERICAL RESULTS 


Flutter calculations have been carried out for the three aeroelastic test 
cases listed in Table 1, and compared to previously published calculations based 
on various TSD codes. Case A is the same as studied by Isogai (Refs. 5,6) and 
later by Edwards et al. (Ref. 7) and also by Weatherill and Ehlers (Ref. 8). 
Note that the elastic axis location "a" is ahead of the leading edge; the idea 
here is to simulate the vibratory behavior (in pitch and plunge) of the stream- 
wise sections near the tip of a swept-back wing. Case B has been studied pre- 
viously by Isogai (Ref. 6) and by Ueda and Dowell (Refs. 9,10). Case C was 
introduced by Ueda and Dowell as an example where nonlinear (amplitude) effects 
were clearly discernible, based on LTRAN2 aerodynamics implemented via the 
describing function method. In all cases, the airfoil is forced for 3-6 cycles 
in pure torsion at a reduced frequency of interest, released, and the aero- 
elastic equations are integrated forward in time for another 3-6 cycles. The 
flutter boundary is located by calculating the logarithmic decrement 6 of the 
transient solutions, and interpolating to 6=0 between adjacent solutions with 
different U/b« a . 


TABLE 1 


Aeroelastic Test Cases 


Case 

A 

B 

c 

Airfoil (s ) 

NACA 64A010 

NACA 64A010 

NACA 64A006 

M 

0 

1 

o 

0.80 

0.86 

00 

a 

-2.0 

-0.3 

CO 

o 

1 

x a 

1.8 

0.5 

0.5 

r 2 

3.48 

0.49 

0.49 

a 


60 

60 

60 

V"a 

1.0 

0.2 

0.2 

Refs . 

5, 6, 7, 8 

6,9 

10 
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MESH GENERATION 


The unsteady Euler calculations are carried out on a C -mesh of quadrila- 
teral elements, generated by means of a square root transformation followed by 
selective stretching to compress the grid near the trailing edge. A near field 
view of the resulting mesh is shown in Fig. 2. In the far field, the mesh 
extends to 15-100 chords, depending on direction. The mesh moves with the air- 
foil as a rigid body, i.e. without deformation. Flutter calculations published 
earlier by the authors (Ref. 4) were carried out on a 96 x 16 C-mesh, which 
was found to give adequate engineering accuracy in most, but not all, of the 
cases studied. In the present study, additional calculations have been performed 
on both 96 x 16 and 192 x 32 C-meshes, and the results of Ref. 4 have been 
updated where appropriate. 



FIGURE 2 
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FLUTTER BOUNDARIES FOR CASE A 


Previous studies (Refs. 58), which have been based on a number of dif 
ferent transonic small disturbance (TSD) theories, generally agree that the 
flutter boundary for Case A exhibits a significant "transonic dip", as shown in 
Fig. 3. Also shown in this figure are the results of flutter calculations based 
on the present Euler code, and using an initial forcing amplitude of 0.1 degree 
in pitch. Overall, the agreement with previous TSD calculations are fairly 
good. However, the Euler calculations appear to shift, the bottom of the 
"bucket" toward higher Mach numbers. It is interesting to note that the bend- 
back of the flutter boundary around M - 0.88 observed by Edwards et al . (Ref. 
6) and Weatherill and Ehlers (Ref. 8), is also predicted by the present Euler 
calculations. Not surprisingly, the precise location of the nose of the curve, 
where the flutter boundary has a vertical tangent, was found to be sensitive to 
the mesh size used in the calculations. 



MACH NUMBER 


• Present code; 96 x 16 mesh 

# Present code; 192 x 32 mesh 
O HYTRAN2 (Ref. 7) 

V 0PTRAN2 (Ref. 8) 

A EXTRAN2 (Ref. 6) 


FIGURE 3 


483 



FLUTTER FREQUENCIES FOR CASE A 


Figure 4 shows the flutter frequencies vs. Mach number, with some com- 
parisons to earlier 0PTRAN2 calculations by Weatherill and Ehlers (Ref. 8). As 
first noted by Isogai (Ref. 5), the flutter mode is essentially the first 
(predominantly bending) natural mode. The flutter frequency is close to the 
first coupled natural frequency Wj/(«) a until the nose of the bend-back is en- 
countered. At this point, the flutter frequency increases to a value between 
the two coupled natural frequencies and the flutter mode also changes, although 
it is still associated with the first predominantly bending branch. 



FIGURE 4 
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TRANSIENT SOLUTIONS FOR CASE A AROUND LOWER FLUTTER BOUNDARY AT N « 0.9 


Typical transient solutions for Case A are shown in Figs. 5,6. At a Mach 
No. of 0.9, multiple flutter solutions occur due to the bend back of the flutter 
boundary (see Fig. 3). Figure 5 illustrates the dynamic behavior of the air- 
foil, plotted as h(t)/b and a(t) vs. time, immediately above and below the 
lower flutter point at M = 0.9 In this case, the airfoil is stable below 
(bottom figure) and unstable above (top figure) the neutral stability boundary 

(CL U„/bw vs. M.) Here, the airfoil has been forced for 3-6 cycles in pure 
F F a 

pitch, with an amplitude of 0.1°, and then released at t=0 . 
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TRANSIENT SOLUTIONS FOR CASE A AROUND UPPER FLUTTER BOUNDARY AT M * 0.9 


In the vicinity of the upper flutter point at M = 0.9, the stability 
behavior is reversed from that observed around the lower point. The airfoil is 
now stable for values of nondimensional airspeed U/b <o a , above the neutral 
stability boundary, as shown in the bottom diagram of Fig. 6. Conversely, the 
airfoil is unstable for values of U/bw a below the flutter boundary. 




TIME 


FIGURE 6 
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FLUTTER CALCULATIONS FOR CASE B 


In Table 2, the results from our present flutter calculations for Case B 
are compared to predictions by previous researchers using various TSD codes. 
This case is the same as Case B considered by Isogai in Ref. 6. Note that the 
present Euler calculations predict a somewhat higher flutter speed than the TSD 
calculation by Isogai, but still below the speed predicted by classical linear 
theory. The flutter speed predicted by Ueda and Dowell (Ref. 9), using the 
describing function method based on LTRAN2 aerodynamics, is significantly below 
the predictions of the Euler code. 


TABLE 2 

Comparison of Predicted Flutter Speed for Case B 


Method 

a 

Up/b(*) a 

2kp 

Present 

0.1° 

3.43 

0.201 

Ueda & Dowell 9 

0.25° (4>i ) 

2.95 

0.221 

Isogai 6 

0.1° 

3.25 

0.215 

Linear Theory 

- 

3 . 86 

0.210 
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NONLINEAR DEPENDENCE OF FLUTTER SPEED 
ON INITIAL FORCING AMPLITUDE FOR CASE C 

In Refs. 9-10, Ueda and Dowell investigated the nonlinear amplitude 
dependence of the flutter boundary for Case C, Table l. They found a distinct 
drop in the flutter speed as the amplitude of the effective induced angle of 

attack <|>i = a + h c /U exceeded about 1°, where h c is the plunging velocity at 
midchord. Figure 7 shows results from the present Euler calculations, plotted 
as flutter speed vs. initial forcing amplitude in pitch (prior to release). 
Note that the flutter boundary is not very sensitive to a in the range 0°-5°, 
and that the results obtained are sensitive to the initial forcing frequency. 




FIGURE 7 
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TRANSIENT SOLUTIONS FOR CASE C 


Typical stable and unstable transient, solutions are shown in Figs. 8 and 9, 
corresponding to initial forcing amplitudes of a = 1° and 4°, respectively. 
The flutter mode is again a predominantly bending mode and emerges quickly 
(within a couple of cycles), despite the fact that the initial disturbance is 
pure torsion. This rapid convergence toward the significant aeroelastic mode 
was also observed in most of the transient solution of Cases A and B as well. 



TIME TIME 


FIGURE 8 


FIGURE 9 
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CONCLUSIONS 


1. Typical section flutter calculations based on the two-dimensional unsteady 
Euler equations are now feasible. 

2. Flutter speeds predicted by the present Euler code are in good overall 
agreement with previous TSD calculations, except in cases where strong 
shocks are present. 

3. The Euler code calculations predict a transonic dip similar to the corre- 
sponding dips predicted by TSD codes, but shifted toward higher Mach numbers. 

4. Multiple flutter points occur at certain Mach numbers, caused by a bend- 
back of the flutter boundary. 

5. The amplitude dependence of Up appears to be less than might be expected. 
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MODEL AND TEST CONDITIONS 


The supercritical airfoil on which the measurement was made, Sc(2)-0714, was developed 
at Langley (Ref. 1). was fourteen percent thick, and had a six-inch chord and an eight inch 
span. The model was machined from Vascomax-200 which has superior dimensional 
stability properties at cryogenic temperatures. The tests were conducted in the Langley 0.3 
Meter Transonic Cryogenic Tunnel (Ref. 2), ,3-m TCT, at Reynolds numbers, R, which 
varied from 6 x 10 6 to 35 x 10 6 at Mach numbers between 0.65 and 0.74. The higher 
Reynolds numbers were near the edge of the tunnel Operating boundary which has a 
stagnation temperature of 120° Kelvin (-243°F) and a stagnation pressure 6.5 atm (Ref. 

2). This tunnel was used in the Advanced Technology Airfoil Test (Ref. 3) program in 
extensive steady flow airfoil studies that demonstrated the necessity for high Reynolds 
number testing. 


EXTERNAL VIEW OF MODEL 



ORIGINAL PAGE 

494 BLACK AND WHITE PHOTOGRAPH 



MODEL CONFIGURATION 


A cavity machined in the underside of the airfoil provided the space necessary to mount the 
transducers. The cavity was closed by a cover plate on which some lower surface 
transducers were mounted. The wing was supported on one end by a close-fitting tang fixed 
to a driving plate with machine screws; this end, on the left of the figure, was sealed with 
epoxy . The other end was supported by an integral shaft which rotated in a bushing in the 
tunnel side wall plate. A sliding seal of felt was used to seal the gap between the end of the 
oscillating airfoil and the fixed tunnel sidewall plate. The position of the supports was 
designed to locate the pitch axis at thirty-five percent chord. 


INTERNAL CONFIGURATION OF MODEL 
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TRANSDUCER LOCATIONS 
(TOP SURFACE) 


Forty-three unsteady pressure transducers were mounted internally in the model. Forty of 
the transducers were mounted in receptacles connected to the orifice by a short length of 
tubing. The remaining three transducers, close mounted, were mounted with the transducer 
head less than 0.1 inch below the surface of the wing. The distribution of the twenty-seven 
upper surface transducers is shown in the figure. The receptacle mounted transducer 
orifices were aligned alternately in two rows 0.25 inches on either side of the airfoil center 
line. The close mounted transducers orifices and reference orifices were located 0.5 inches 
from the center line. The orifices of the close-mounted transducers were paired with 
receptacle-mounted transducer orifices for comparison purposes. The orifices were 
distributed every 2% chord to x/c of 0.1 and 4% chord to x/c of 0.7. 
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TRANSDUCER LOCATIONS 
(BOTTOM SURFACE) 


This figure shows the distribution of the 16 transducer orifices on the lower surface; the 
orifice at the leading edge measured only static pressure. The distribution of the 1 5 
receptacle mounted transducers is more sparse than on the upper surface and is 
concentrated in regions of largest pressure gradient and is 2% to an x/c of 0.1 and 0.5 
thereafter. 
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ELEMENTS OF THE TRANSDUCER SYSTEM 


The system consisted of transducers, designed for cryogenic application, with a 10 psi range 
and with outputs of between 5 and 9 mv/psi. Each transducer was mounted in a receptacle 
which in turn was connected to the 0.015 inch diameter orifice by a 0.75 inch length of 
0.030 inch i.d. tubing. Each transducer was referenced to a manifold which in turn was 
vented to one of five static reference orifices. A reference transducer measured the 
pressure differential between the manifold and the tunnel static pressure. The connection 
between the manifold and the reference orifice was interrupted by a porous flow restrictor 
which damped out the oscillating pressure from the reference orifice. A series of tests were 
conducted before the model was fabricated to examine the effects of orifice diameter, tube 
diameter and tube length on the dynamic response of the system. At atmospheric conditions 
there was no significant reduction of dynamic amplitude response or phase shift of the test 
configuration up to 100 Hz. 



Restrictor case 
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MODEL INSTALLATION M 


The large variations in temperature (120° k to 320° K) and stagnation pressure (1.4 atm. 
to 6 atm.) over the operating range of the 0.3-m TCT results in plenum wall deformations 
that required special consideration in the design of the oscillating drive system. The 
photograph of the test section, with the ceiling removed, shows the installation of the' airfoil 
and drive system. The hydraulic-rotary actuator required the maintenance of precise 
alignment during the test. Since the test section floats on a cable suspension system to 
accommodate thermal contraction, the actuator and supporting structure were also supported 
by a system of cables and counterweights to enable them to move with the plenum wall. 
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SCHEMATIC OF MODEL INSTALLATION 


The critical elements of the oscillating system are identified in this figure. The hollow 
aluminum drive shaft had fixed point supports at the rotary shaft and at a Teflon bushing and 
a pressure seal located on the tunnel plenum wall plate. The shaft was attached to the 
rotating sidewall wall drive disk through a bellows that allowed in-line shaft movement with 
the plenum wall. The rotating drive disk was Teflon coated on its circumferential bearing 
surfaces and had a slot to accommodate the wing tang. The tang was hollow to provide a path 
for transducer cable and tubing. The other edge of the wing was supported by an integral 
hollow shaft and a bushing in the sidewall plate. The hollow shaft allowed a path for the 
remaining transducer cables. The rotary transducer, attached to the shaft, was heated with 
surface heaters under thermostat control and the assembly was covered by an insulating can. 
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DATA ACQUISITION AND REDUCTION 


Steady Pressures - The airfoil and tunnel instrumentation signals were fed to the tunnel 
data acquisition system through a 10 Hz low-pass filter, digitized at 20 samples/sec and 
averaged over a one second interval. 

Unsteady Pressures - The signals from the amplifiers were recorded on two 28 channel 
recorders operating at 15 inches per second. This analog data was then digitized at 32 
samples/cyc of oscillatory motion for 64 cycles and the harmonic components of the 
unsteady pressure were determined from FFT analysis. These components were normalized by 
the harmonic pitch amplitude in degrees. All phase angles were relative to wing position. 

Tunnel Corrections - Sidewall boundary-layer and angle of attack corrections were applied 
to the measured steady pressure results. The sidewall boundary-layer corrections are based 
on the analysis of Ref. 4 which is used in Ref. 5 with measured values of sidewall 
displacement and momentum thickness to compile the tables which were used to correct the 
experimental values in this paper. The angle of attack corrections described in Ref. 6 
(sometimes referred to as the "Barnwell-Davis-Moore" correction) adjust the analysis of 
Davis-Moore with experimental data. The wall induced downwash over the airfoil in the 
0.3-m TCT for Cl = 1.0 is: 


5a = 1.73245° 


i 
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MACH NUMBER AND REYNOLDS NUMBER TEST CONDITIONS 


The test was designed to explore the effects of Reynolds number on unsteady pressures and to 
generate a data base for validating unsteady-aerodynamic computer codes. The test 
conditions as defined by Mach number and Reynolds number are shown in the figure. Test 
points were taken at the design Mach number of 0.72, determined from data of a previous 
test, at test Reynolds numbers varying from 6 x 106 to 35 x 10 6 . A total of 976 test 
points were taken. The primary data base was tak$n for pitch-oscillation frequency between 
5 Hz and 40 Hz at an amplitude of ±0.25° as indicated by the open and solid symbols. Once 
this data was in hand, the pitch amplitude was increased to ±0.5° and ±1 .0° and the pitch 
frequency increased to 60 Hz at test conditions indicated by the solid symbols. 
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COMPARISONS OF STEADY TEST RESULTS WITH CALCULATED 
RESULTS AT A TUNNEL MACH NUMBER OF 0.72 


a t « 2.5 degrees 


The next four figures give the steady pressure distributions for four angles of attack at two 
Reynolds numbers, 6 x 10 6 and 30 x 10 6 . Experimental data, shown as symbols, are compared 
with calculated results from GRUMFOIL computer code (Ref. 7) which are shown as solid lines. 
The GRUMFOIL code consists of a full-potential equation flow solver integrated with a viscous 
boundary layer model and may be entered by specifying either a or C|. The corrected values of 
Mach and C| were used as the input data for the computed results. Below each figure are 
listed M, a , and C| for the tunnel test conditions, the corrected values, and the values resulting 
from the GRUMFOIL calculations. 


a j = 2.5 deg 

o Upper surface o Upper surface 



x/c x/c 

Tunnel Corrected Grumfoil Tunnel Corrected Grumfoil 

M 0.720 0.701 0.701 M 0.719 0.705 0.705 

o 2.504 0.844 1.385 a 2.51 0.756 1.399 

C 1 0.9581 0.9753 0.9837 C ] 1.0123 1.0256 1.0336 
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o pi 


COMPARISONS OF STEADY TEST RESULTS WITH CALCULATED 
RESULTS AT A TUNNEL MACH NUMBER OF 0.72 

a { - 2.0 degrees 

The comparisons between the experiment, shown as symbols, and GRUMFOIL calculations, solid 
lines, are very good. The shock moves aft by approximately, 8% to 1 0% of chord for a given 

value of tunnel mean angle of attack, a t . when Reynolds number is increased from 6 x 10 6 to 
30 x 106. The code under-predicts the position of the shock at both Reynolds numbers by 
approximately 2-3% of chord even through Cl is matched. 


COMPARISONS OF STEADY TEST RESULTS WITH CALCULATED 
RESULTS AT A TUNNEL MACH NUMBER OF 0.72 

a t =2.0 deg 

° Upper surface o Upper surface 



0 .2 .4 .6 .8 1.0 0 .2 .4 .6 .8 1.0 


x/c x/c 

Tunnel Corrected Grumfoil Tunnel Corrected Grumfoil 

M 0.72 0.701 0.701 M 0.721 0.705 0.705 

2.002 0.525 0.961 a 1.997 0.393 0.910 

0.8523 0.8676 0.8757 C ] 0.926 0.939 0.9453 
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COMPARISONS OF STEADY TEST RESULTS WITH CALCULATED 
RESULTS AT A TUNNEL MACH NUMBER OF 0.72 

Qj ~ 1.5 deg 



x/c x/c 

Tunnel Corrected Grumfoil Tunnel Corrected Grumfoil 

M 0.719 0.701 0.701 M 0.718 0.705 0.705 

a 1.495 0.201 0.493 a 1.501 0.051 0.552 

C 1 0.7467 0.7601 0.7680 C 1 .8373 0.849 0.8555 


COMPARISONS OF STEADY TEST RESULTS WITH CALCULATED 
RESULTS AT A TUNNEL MACH NUMBER OF 0.72 

a t ~ 0.0 deg 



x/c x/c 


M 

Tunnel 

0.720 

Corrected 

0.701 

Grumf oi 1 
0.701 

M 

Tunnel 

.721 

Corrected 

0.705 

Grumf oi 1 
0.705 

a 

0.004 

-0. 92 

-0.398 

a 

-0.005 

-1.036 

-0.715 

C 1 

0.5288 

0.5383 

0.5484 

C 1 

0.5951 

0.6034 

0.6099 
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COMPARISON OF LIFT COEFFICIENT VERSUS CORRECTED 

ANGLE OF ATTACK 


Lift coefficients for several cases are shown in this figure plotted against corrected angle of 
attack and against angle of attack as computed by GRUMFOIL code for input values of Mach 
number and C| for Reynolds numbers of 6 x 106 and 30 x 106 . The angles calculated by the 
code are consistently larger than those determined by the correction procedure of Ref. 6. 
Irrespective of angle of attack corrections, an increase in Ci of approximately 0.1 is shown 

as Reynolds number is increased from 6 x 106 to 30 x 106 . This increase results from the 
rearward movement of the shock shown in the previous figures. 
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UNSTEADY PRESSURE TEST RESULTS AT A TUNNEL 
MACH NUMBER OF 0.72 AND AT a = 0.25 DEGREES 


a t = 2.05 degrees 


The effect of Reynolds number and frequency of pitch oscillation on the upper surface unsteady 
pressure distribution is shown in the next two figures. Results are given in terms of the 
modulus of the unsteady pressure coefficient normalized by the oscillating pitch angle, a, and 
the phase angle, <|>, between the unsteady pressure and the oscillating pitch angle. Results are 
shown for a m 1 and 2° at R - 6 x 106 and 30 x 10 6 for two oscillation frequencies, 5 Hz and 
20 Hz at a pfrch amplitude of ±0.25°. 
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UNSTEADY PRESSURE TEST RESULTS AT A TUNNEL 
MACH NUMBER OF 0.72 AND AT at = ±0.25 DEGREES 

a ( = 1.04 degrees 


The shock wave, identified by the peak in the unsteady pressures, moves aft about 8% to 
10% chord as R is increased from 6 x 10 6 to 30 x 10 6 at the same tunnel test jngles. The 
unsteady pressures, at both R, are significantly greater ahead of the shock ata { = 1° than 

at 2o, but there is no significant difference in the pressure modulus due to the change in 
frequency from 5 Hz to 20 Hz. For both angles of attack and R the pressures ahead of the 
shock are approximately 180° out of phase with the wing oscillation. At the shock the phase 
angle abruptly changes from -180o to 0°. Behind the shock the phase angle is erratic at a t 
= 1 .oo and is more dependent on frequency than at a = 2.0o. 


UNSTEADY PRESSURE TEST RESULTS AT A TUNNEL 
MACH NUMBER OF 0.72 AND ATa=±0.25 
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VARIATION OF |Cp| /DEGREES AND <f> WITH PITCH AMPLITUDE 
AT M = 0.72 AND R = 30 X 106 
f = 40 HZ 


The effect of varying the amplitude of oscillation at M = 0.72 and R = 30 x 106 is shown in the 
next two figures. Upper surface pressure modulus and phase are shown for three pitch 
amplitudes (0.25, 0.5, 1.0°) at mean angles of 1.0 and 2.0° at f = 40 and 60 Hz. A reduction 
and broadening of the shock-generated peak amplitude is evident as the pitch amplitude, a, is 
increased at both frequencies and mean angles. 
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VARIATION OF |Cp| /DEGREES AND <t> WITH PITCH AMPLITUDE 
AT M s 0.72 AND R = 30 X 106 
f = 60 HZ 


A secondary peak in the modulus of the oscillating pressure is evident immediately behind the 
shock in this and the other unsteady pressure distribution figures. The amplitude of this second 
peak is greater at a =1.0° than at a ( = 2.0°. The phase angle between pressure and wing 
motion is approximately -1 80° between the leading edge and immediately behind the shock, at 
which point there is a sudden decrease to below -80°. There is less deviation in phase 
for the “ = 2.0° data than for the ” =1 .0° data. 
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SONIC REGIONS AT M = 0.72 CALCULATED BY GRUMFOIL CODE 


The secondary peak in the pressure modulus could be attributed to flow separation and 
reattachment as discussed in Ref. 9. However, an inviscid calculation using XTRAN2L 
(Refs. 10,11) computer code predicts this secondary peak in the same relative location. The 
sonic regions calculated from GRUMFOIL suggests that a more probable reason for the 
secondary response derives from the supersonic region above the airfoil. At a = 1 -0° 
there is a secondary supersonic region behind the shock which is engulfed by the primary 
supersonic region when the angle of attack is increased to a ( = 2.0°. Tijdeman (12) and 
others have noted that the flow in the supersonic region prior to the formation of a shock is 
characterized by a substantial increase in unsteady pressure. 
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VARIATION OF |Cp| /DEGREES AND $ WITH FREQUENCY 
AT M = 0.72 AND R = 30 X 106 
a = 0.25 DEGREES 


The effect of frequency on the modulus and phase of the upper surface unsteady pressures is 
shown in the next two figures for M = 0.72 and R = 30 x 10 6 . The measurements are shown 
for frequencies of 5, 15, 40, and 60 Hz at mean pitch angles of 1 and 2 ° and at amplitudes of 
0.25 and 0.5°. In general the excursion of the shock is reduced at 60 Hz and again the 
amplitude of the second peak is reduced for all frequencies as a is increased from 1 .0 to 
2.0o. 


ICpI /deg 
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VARIATION OF |Cp| /DEGREES AND $ WITH FREQUENCY 
AT M = 0.72 AND R = 30 X 106 
a = 0.5 DEGREES 

The phase angle shows a dependency on frequency in detailed sense having the same overall 
characteristics as in the previous figures. The phase angle between the pressure and the 
airfoil motion is approximately - 180° from the leading edge to immediately behind the 
shock where it increases rapidly to approximately 0°. 



x/c 


513 



TIME HISTORIES AT FIVE CHORD STATIONS FOR 
A = 0, M = 0.72, AND R = 35 X 106 

a = 0 degrees 

The pressure transducers used to measure unsteady pressures generated by airfoil pitch 
oscillation were also used to measure unsteady pressures when the wing was held at a fixed 
angle of attack. The next two figures show time histories at five chord stations whose 
locations are shown by solid symbols on the plot of the static pressure distribution, on the 
right of the figures, for the angles of attack (0 and 2°) being considered. The time histories 
shown in these figures were all taken at a gain of 10, but the transducer sensitivity, given 
with each trace, has not been applied to put the data in engineering units. 


A, x/c =0.14, Sens. -8.45 mv/PSI 

B, x/c = 0.26, Sens.- 7.23 mv/PSI 

llUyJUML^JU !■!»' 

C, x/c=0.46, S ens. -8. 4 mv/PSI 


D, x/c = 0.62, Sens. = 5.83 mv/PSI 



E, x/c =0.75, Sens. = 6.85 mv/PSI 


i i i i 1 1 1 1 



0 .2 .4 .6 


.8 1.0 1.2 1.4 



Time, sec 
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TIME HISTORIES AT FIVE CHORD STATIONS FOR 
A = 0, M = 0.72, AND R = 35 X 106 

a = 2 degrees 

At a ( = 0° the time histories have the characteristics of a turbulent boundary layer and 
are in sharp contrast with the time histories of a t = 2°. At a t of 2° the pressure is 
quiescent at x/c of 0.14 and 0.28, at x/c of 0.46 the effect of shock movement on the 
pressure response is observed, which increases at the foot of the shock, x/c = 0.62, where 
turbulence is also apparent. The observable differences between the flows at the two angles 
of attack are the more favorable pressure gradient and the presence of a shock wave at a = 

2°. The time histories indicate that laminar flow was present at a t = 2° and that 
transition to turbulence was between an x/c of 0.28 and 0.46 corresponding to transition 
Reynolds numbers between 9.8 x 10 6 and 16.1 x 10 6 . The possibility exists that long runs 
of laminar flow existed intermittently during the tests. 


A, x/c = 0.14. Sens. = 8. 45 mv/PSI 



B, x/c = 0.26, Sens. = 7. 23 mv/PSI 



0 .2 .4 .6 .8 1.0 1.2 1.4 


Time, sec 
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CONCLUSIONS 


Steady and unsteady pressures have been measured on a 14 percent supercritical airfoil at 
transonic Mach numbers at Reynolds numbers from 6 x 10 6 to 35 x 106 . Instrumentation 
techniques were developed to measure unsteady pressures in a cryogenic tunnel at flight 
Reynolds numbers. Experimental steady data, corrected for wall effects show very good 
agreement with calculations from a full potential code with an interacted boundary layer. 
The steady and unsteady pressures both show a shock position that is dependent on Reynolds 
number. For a supercritical pressure distribution at a chord Reynolds number of 35 x 10 6 
laminar flow was observed between the leading edge and the shock wave at 45 percent chord. 
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MODEL WING IN WIND TUNNEL 


The model wing consists of a set of fiberglass panels mounted on a steel spar that 
spans the 8 ft. test section of the UTRC Large Subsonic Wind Tunnel. The pitch angle 
of the model is set by hydraulic actuators attached to each end of the spar and 
controlled using a dual channel closed-loop system. This allows 1024 samples of each 
data channel to be acquired at identical positions during each cycle of the wing 
motion. The first use of this system was to measure surface pressures and flow 
conditions for a series of constant pitch rate ramps and sinusoidal oscillations at 
Mach numbers between 0.2 and 0.4, Reynolds numbers between 2 and 4x10^, and pitch 
rates between A = cd/2U = 0.001 and 0.02. 
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AIRFOIL SECTION 


A Sikorsky SSC-A09 airfoil with a 17.3 in. chord was used. This is a 9 % thick 
supercritical section designed for low drag at high subsonic Mach numbers. 
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INSTRUMENTATION 


The wing-mounted instrumentation consists of 72 miniature pressure transducers and 8 
surface hot film gages. The pressure transducers are mounted inside the fiberglass 
skin and are connected to the surface by short pipettes. This technique provides a 
point measurement on a smooth surface contour. The primary chordwise pressure 
transducers are located in arrays of 18 on each surface, 0.5 chordlengths from the 
tunnel centerline, and 2.3 chordlengths from the side wall. The transducers are 
arranged in a segmented Gaussian array from x/c = 0.005 to x/c = 0.99. The other 36 
transducers are located in additional arrays at 1.8 and 1.4 chordlengths from the 
wall, intended for use in future swept wing and finite tip experiments. The hot film 
gages are located in an array parallel to the primary pressure array, from x/c = 

0.026 to x/c = 0.88. Detailed calibrations over both pressure and temperature ranges 
make it possible to measure both mean and unsteady pressures to within 1 % of the 
calibration range. An interactive data system is used to acquire single samples of 
time histories, or to obtain an ensemble average of the time histories, based on 20 
cycles of the motion. These data are digitized, converted to pressure coefficient, 
and stored on magnetic tape. 
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TEST CONDITIONS 


Unsteady data were acquired for 36 constant pitch rate ramps and 9 sinusoidal 
oscillations, as shown in the table. The maximum pitch rate of 360 deg/sec occurred 
at A * 0.02 and M * 0.2, and the minimum of 18 deg/sec at A * 0.001. The maximum 
rate is lower than the maximum reached in many smaller scale experiments (Refs. 

1-5), but is larger than that for both the "typical" (A = 0.001) and "minimum time" 
(A - 0.0044) maneuvers described in Ref. 6. The wing angle of attack was limited by 
the support system to a maximum of 30 deg, which was less than the maximum obtained 
in the smaller scale tests, but sufficient to include all of the primary 
stall-related events at the pitch rates used. 


STEADY 

M * 0.2, -5°<a<28° 

M = 0.4, 0°< a <20° 

UNSTEADY RAMPS 



UNSTEADY OSCILLATIONS 


k = Cu/2U 


M 

a 

.025 

.050 

.100 

0.2 

5° + 5°sinut 

X 




10° + 10°sinuit 

X 

X 

X 


20° + 10° si nut 

X 

X 

X 

0.3 

9° + 8°sinut 


X 



12° + 8°sinu)t 



X 


(9 SINUSOIDS) 
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STEADY PRESSURES 


The steady pressure distributions have a strong leading edge suction peak prior to 
stall. Maximum suction occurs at a = 13 deg and reaches values of Cp = -8.2 at 
M * 0.2, and C = -5.2 at M =0.4. The sonic pressure coefficient at M = 0.4 is 
-3.7, so a small supersonic zone exists at the leading edge, terminated in a shock 
near x/c - 0.03. The compressibility effects associated with this zone create the 
lift curve differences seen earlier. 



0.0 0.2 0.4 0.6 0.8 1.0 

X/C 
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RAMP TIME HISTORIES 


Airfoil motion for a ramp consists of an initial delay 
minimum angle, a constant rate increase to the maximum 
maximum angle, and a return to the initial condition, 
a small portion, T, of this cycle. Time histories of 
acquisition period are shown in the figure for several 
The ramp begins at nondimensional time, r = t/T, of 0. 
pitch increase is quite linear with time and has sharp 
higher pitch rates the damping of the hydraulic system 
maintains a nearly linear pitch rate near stall. 


of several seconds at the 
angle, a second delay at the 
Data are acquired only during 
the pitch angle during the data 
0 to 30 deg ramps at M = 0.2. 
125 and ends at t = 0.625. The 
corners for A < 0.005. At 
rounds the corners, but still 


= 0.0->30.0 deg 

= .20 

A= 0.0010 
A* 0.0025 
A= 0.0050 
A= 0.0100 
A= 0.0200 


0.0 

0.2 

0.4 

0.6 

0.8 

1.0 



Time 
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UNSTEADY PRESSURES, UPPER SURFACE 


Ensemble averaged upper surface pressure time histories for 11 of the 18 measurement 
stations are shown for a 0 to 30 deg constant pitch rate ramp at A = 0.005, M = 0.2. 
The ordinate scale on the lower left refers absolutely to the x/c - 0.005 curve, and 
all other curves are vertically offset by AC - 1.0. Each of these curves is 
referenced to its own origin (upper left scale). The pressure responds smoothly to 
the imposed pitch angle until maximum suction is reached at r = 0.4 (a = 15.5 deg) at 
letter A. Peak suction at x/c - 0.005 is Cp - -12.5, corresponding to a local Mach 

number of 0.84. The rapid local increase in suction associated with passage of the 
vortex (letter B) is only observed for x/c > 0.1,, implying that the vortex forms 
ahead of this position and travels downstream thereafter. The vortex pressure 
propagates along the chord at 0.16U, and reaches the trailing edge at r * 0.45. For 
t> 0.50 the upper surface separation is massive and the pressures are virtually 
invariant in both time and position (letter C). All events are completed well before 
the ramp ends at r = 0.625. 
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UNSTEADY PRESSURES, LOWER SURFACE 


The lower surface does not separate, and responds primarily to the migration of the 
stagnation point near the leading edge (from motion initiation at t = 0.125 to r* 
0.42), and to the upper surface massive stall near the trailing edge (letter A), 
where the pressure is required by continuity to match that on the upper surface. In 
addition, a periodic oscillation in these ensemble averaged pressures is detected 
near the lower surface trailing edge (letter B). The oscillation frequency of 62Hz 
is very close to the 65Hz frequency calculated for the vonKarman vortex street shed 
by a circular cylinder having a diameter, D, equal to the vertical projection of the 
airfoil chord, c sin a. From Ref. 7, a cylinder will generate a vortex street at a 
Strouhal number, Df/U, of 0.21 for Reynolds Numbers less than 10^ and at fD/U - 0.27 
for Re > 3 x 10^. No regular vortex street is formed when the Reynolds number is 
between these limits. In the present case the Reynolds number of the cylinder is 
equivalent in projected area to the airfoil at M - 0.2 and a - 30 deg is 10^, near 
the boundary where periodic oscillations should cease. 
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INDIVIDUAL RECORDS, LOWER SURFACE 


A larger post-stall oscillation amplitude is found for individual records of the 
motion than for the ensemble average. The top time history in this figure is the 
ensemble averaged trailing edge plot from the previous figure, and the remaining 
curves are a set of (nonconsecutive) individual records at the same location and for 
the same set of flow conditions. Each curve is offset vertically by Cp — 0.5. The 
oscillations have the same frequency and a similar maximum amplitude during each 
record, but they are not well correlated in phase. The resulting cancellation 
reduces the ensemble averaged amplitude. A separate oscillation at a lower 
frequency is often present immediately following the passage of the stall vortex at 
t - 0.45 (letter A), and may result from secondary stall vortices. 
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HOT FILM TIME HISTORIES 


Hot film gages provide additional information on the surface flow conditions, as 
shown by the six time histories in this figure. The data are self-scaled, so each 
varies over one unit. As indicated on the leading edge trace, the boundary layer is 
initially laminar at the first two stations, goes through transition, becomes 
turbulent, and eventually becomes separated. From x/c = 0.192 and aft the flow is 
always either turbulent, or separated, and is never laminar. Ahead of the stall 
vortex release point (x/c * 0.1) the heat transfer drops when the boundary layer 
separates at r = 0.4. Downstream of the release point the heat transfer drops 
slightly at separation but then rises rapidly as the high velocities induced by the 
stall vortex pass each gage. There is excellent agreement between the times of 
separation as measured by the hot film gages and by the adjacent pressure 
transducers. 
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UNSTEADY PRESSURE, UPPER SURFACE, MAX PITCH RATE 


At the maximum pitch rate of A = 0.02 the upper surface pressure time histories are 
qualitatively similar to the lower pitch rate results prior to stall, but are 
significantly different after stall. The unsteady increments to the airloads are 
strongly influenced by the strength of the stall vortex, which may be estimated by 
measuring the local change in Cp at some point along the chord downstream of the 
vortex release point, say at x/c - 0.302. This increment increases from AC - 1.0 
at A - 0.001 (not shown here) to AC - 1.5 at A - 0.005 and to AC - 2.5 at P A - 
0.020 (below). p p 
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VORTEX PROPAGATION TRACES 


'The vortex propagation speed along the airfoil may be estimated using the times of 
minimum pressure at each transducer. These points are plotted below for 5 pitch 
rates. The region of constant vortex speed generally begins near x/c = 0.10 and ends 
between x/c = 0.80 and 0.90, and the speed in this region increases approximately 
linearly with pitch rates, from 0.13U at A = 0.001 to 0.33U at A = 0.020. These 
results are consistent with those of previous investigations (Ref. 8 for sinusoidal 
motions) . 
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UNSTEADY FORCE AND MOMENT COEFFICIENTS 


Integrated lift, pressure drag, and pitching moment (parts a, b, c) exhibit 
characteristic behavior with increasing pitch rate: (1) Before the stall vortex is 
formed, increasing pitch rate decreases the lift slope, decreases pitching moment, 
and increases drag. The lift slope effect agrees with the results of Ref. 2 at lower 
Reynolds number. (2) There is a rapid buildup of lift as the leading edge vortex 
forms. The unsteady increment added to the quasi-steady lift increases from AC^ * 

0.4 at A = 0.001 to 1.1 at A = 0.020. (2) The pressure drag increases smoothly 

before stall, rises and falls rapidly as the stall vortex travels over the chord, 
then increases slowly as the pitching of the airfoil rotates the aerodynamic force 
vector. (4) The pitching moment in attached flow becomes more negative at higher 
pitch rate, following the prediction of thin airfoil theory (cf. Eq. 4-171 in Ref. 

9). (5) The peak negative C M after stall increases with pitch rate, from C M = 

0.22 at A = 0.001 to -0.52 at A > 0.020. 
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SEQUENCE OF STALL EVENTS 


The sequence of stall events is systematic and each appears to be linear with pitch 
rate. First the stall vortex forms and and start to rise. Moment stall 
occurs when the vortex is released. Maximum occurs as the vortex travels 
downstream along the chord, and when the vortex reaches the trailing edge, the 
minimum and maximum C D are obtained. At these pitch rates the angle when moment 
stall occurs increases approximately linearly with A, in agreement with results for 
a sinusoid at similar peak pitch rates (Ref. 8). This does not agree with the 
square root correlation postulated by Gormont (Ref. 10). 
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UNSTEADY PRESSURES, UPPER SURFACE, MAX MACH NO. 


A supersonic zone at the leading edge has a profound effect on chordwise pressures at 

freestream Mach number M = 0.4. The time histories for this 0 to 20 deg ramp at A - 

0.005 are not shifted vertically and the ordinate is an absolute measure for all 
curves. At x/c = 0.005 the flow becomes supersonic at r = 0.34 (a = 7.6 deg) and 
reaches a peak local Mach number of 1.3 at r = 0.45 (a* 10.6 deg). The sharp rise 
and fall of the pressure at x/c = 0.026 (letter A) is caused by the movement of the 
shock downstream past this station at r = 0.40 and back upstream at r = 0.46. A 
distinct stall vortex is released at r = 0.46 and x/c = 0.060, just downstream of the 
shock (letter B). Several of the vortex characteristics differ from those at M • 

0.2: (1) the release point is at x/c = 0.060 instead of 0.10; (2) the vortex speed 

is approximately 10% less than the speed at M = 0.2; (3) the strength of the vortex 
is reduced by approximately 50%; and (4) the pressure signature of the vortex is not 
observed downstream of x/c = 0.57. In addition, no clearly defined oscillations are 
present after stall at M - 0.4. The Reynolds number of the equivalent bluff body is 

1.4 x 10^ for this case, above the maximum of 1 x 10^ for a stable vortex street 

(Ref. 7). 
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Lift Force, 







UNSTEADY PRESSURES, UPPER SURFACE, SINUSOID 


ORIGINAL PAGE IS 
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Sinusoidal motion differs from ramp motion in two important aspects: in a sinusoid 

the airfoil does not start from a steady-state condition, and the pitch rate is 
constantly changing. Nonetheless, there are many qualitative similarities between 
the responses to the two motions. For the upper surface time histories shown here 
(a = 20 - 10 coso»t at M = 0.2 and k « 0.05) the motion begins at r = 0.0, reaches 
maximum pitch angle at t « 0.50, and returns to the minimum value at r * 1.0. Each 
time history here is offset vertically by AC = 1.0. The characteristics are 
similar to the constant pitch rate results sRown earlier for A = 0.005, which is 
comparable to the instantaneous sinusoidal pitch rate at stall of A = 0.008. 
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UNSTEADY PRESSURES, LOWER SURFACE, SINUSOID 


The lower surface pressures for this same sinusoid are remarkable in that the maximum 
pressure change occurs at the trailing edge (to match the stalled suction surface 
pressure level), exhibits an initially damped behavior, and then stablilizes into a 
coherent vortex street response. The latter is noteworthy because in this 
ensemble-averaged data set the vortex street frequency must be synchronized with the 
fundamental frequency of the sinusoid. No single records were taken to permit close 
examination of the vortex street oscillation. 
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TRANSITION POINT DETECTION 


Accurate determination of the point of boundary- layer transition is critical to many 
numerical computations. Surface heat transfer gages are usually required to 
determine the unsteady motion of the transition point. These gages are generally 
cumbersome to install and each requires a dedicated active anemometer circuit to 
operate. The techniques developed during the course of this study show that it may 
be possible to locate the transition point using the local increase in rms pressure 
measured by the unsteady pressure transducer as the local flow changes its character. 
An example of the correlation between the hot film and pressure transducer results 
is shown below. All results are at x/c « 0.026 for a 0 to 30 deg ramp at A - 0.02 
and M = 0.20. The vertical pressure scale has been severely foreshortened, and only 
a small portion of the time scale is displayed. The increase in pressure rms 
corresponds to the increase in ensemble averaged heat transfer (upper curve) and 
to the spike in the rms heat transfer (second curve). No change at transition can 
be discerned in the ensemble-averaged pressure (bottom curve). If these promising 
early results are confirmed by additional correlations, this method may make it 
easier to locate transition in complex three-dimensional and unsteady flows. 
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NOMENCLATURE 


A pitch rate, ac/2U, rad/sec 

c airfoil chord, m 

Cp section pressure drag coefficient, D/Qc 

C. section lift coefficient, L/Qc 

L 2 
section pitching moment coefficient about x/c * 0.25, m/Qc 

Cp pressure coefficient, (P-Pstatic)/Q 

k reduced frequency of sinusoidal motion, o»c/2U 

2 

Q freestream dynamic pressure, 0.5 p U , Pa 

t time from start of data acquisition, sec 

T data acquisition period, sec 

Re Reynolds number, c U/i> 

U freestream velocity, m/sec 

x distance from airfoil leading edge, m 

a geometric pitch angle, deg 

2 

v kinematic viscosity, m /sec 

3 

p air density, kg/m 

r nondimensional time, t/T 

a> frequency of sinusoidal oscillation, 2 n f 
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CONCLUSIONS 


Conclusions - Increased Pitch Rate 

• Stall events are delayed 
•Stall vortex is strengthened 
•Vortex propagation speed increases 
•Unsteady airloads increase 


Conclusions - M = 0.3, 0.4 

• Supersonic zone near leading edge 

• Stall vortex is weaker 

• Unsteady airloads are reduced 

Additional Conclusion 

• Post stall vortex shedding when Re bi u ff < 10® 

• Shedding synchronized when A > 0.01 

• Sinusoids and ramps qualitatively similar 


Note: This work was supported by the U.S. Air Force Office 

of Scientific Research under Contract F49620-84-0082 
and will be reported in Ref. 11. 
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PICTURE OF MODEL IN WIND TUNNEL 

A joint Langley-Lockheed wind tunnel test was undertaken involving 
this model. The motivation for this test is explained in the figures to 
come. 
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CONVENTIONAL VERSUS SUPERCRITICAL AIRFOILS - GEOMETRIC SHAPE 


In 1981 Lockheed conducted a wind tunnel test that compared conventional 
and supercritical airfoils while holding stiffness, mass, and planform 
geometric shape constant. 

Model Airfoil Profiles for Instrumented 
Pressure Sections 



CONVENTIONAL AIRFOIL SUPERCRITICAL AIRFOIL 
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CONVENTIONAL VERSUS SUPERCRITICAL AIRFOILS - FLUTTER BOUNDARIES 

This test and other tests show that changing from a conventional airfoil 
shape to a supercritical airfoil shape can greatly reduce the wing's flutter 
speed. The test also showed that there was a region of low damping within 
the flight envelope of this wing with supercritical airfoils. This low 
damping region is shaded in the figure. 
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PROBLEM 

The aerodynamic programs used in flutter analyses do not accurately 
predict the complex flow around supercritical airfoils in the transonic 
flow region. This causes the use of long costly wind tunnel tests and 
empirical weighting factors to modify the analytically predicted flutter 
speeds. The result can be a stiffer, heavier wing than is needed. 

Unsteady transonic aerodynamic programs using Computational Fluid 
Dynamics (CFD) methods show promise of more accurately predicting transonic 
flow, but these programs need to be validated before they can be incorporated 
into a production flutter method. 

To validate the programs, analytical predictions must be correlated 
with steady and unsteady experimental flow data on a flexible, three-dimen- 
sional wing. Most of the data available for correlation is from tests on 
two-dimensional or three-dimensional rigid wings. 

In April 1984, Lockheed-Georgia and NASA-Langley conducted a wind tunnel 
test to obtain all of the types of data needed for CFD program correlation. 
This included steady state data, forced oscillation data, and oscillatory 
data during flutter. 

*r Mew Technologies 

I— lave l_owen F 1 u t t e n Speeds 

*4 P n e s e n t Analytical Methods 

Ane Mot Accunate 

->4 Computational Fluid Dynamics (CFD) 

>4 V e n i f y CFD Pnognams 

*4 Lack: of Test Data fon Connelations 
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TEST OBJECTIVES 


There were three main objectives for this test. 

a) obtain aerodynamic data during flutter for CFD program correla- 
tion 

b) obtain a better understanding of supercritical wing flutter 

c) evaluate the effects that pylons and engines have on wing 
unsteady aerodynamics 

• OBTAIN CONGRUENT FLUTTER AND AERO DATA FOR 
ANALYSIS CORRELATION 


• OBTAIN BETTER UNDERSTANDING OF SUPERCRITICAL 
WING FLUTTER 


• EVALUATE EFFECTS OF PYLONS AND ENGINES ON WING 
AERO DATA 


MODEL CONFIGURATIONS 


Four different model configurations were tested. 

a) stiffer spar, bare wing 

b) nominal stiffness spar, bare wing 

c) nominal stiffness spar with mass simulated engines 

d) nominal stiffness spar with aerodynamic simulated engines 

The first configuration was used only for obtaining forced response 
oscillatory data. The other three configurations were used for obtaining 
both forced response oscillatory data and oscillatory data during flutter. 

SUFFER WING (FOUR TIMES NOMINAL STIFFNESS) 

• BARE WING 

NOMINAL STIFFNESS WING 

• BARE WING 

• WING PLUS DUMMY NACELLES AND PYLONS 

• WING PLUS DUCTED NACELLES AND PYLONS 
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WING PLANFORM AND INSTRUMENTATION LAYOUT 


The wing had an aspect ratio of 7.84. It was constructed on a single 
aluminum spar with a supercritical airfoil. It had eleven mass ballasted 
sections. Five bending and torsion strain gage bridges and five pairs of 
accelerometers were distributed along the wing’s span to define the steady 
and unsteady position of the wing. Instrumentation sections were located 
at 49.6 % and 82.1% span. Each instrumentation section contained 17 delta 
pressure transducers and 7 upper surface pressure transducers from the lead- 
ing to the trailing edge. 

Aeroelastic Model Wing Planform 
and Instrumentation Layout 
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PICTURE OF INSTRUMENTATION SECTION 

This picture shows an opened up instrumentation section. The holes 
for the pressure transducers are visible on the wing's surface. Wires from 
the pressure transducers and the wing's spar are visible inside the model. 
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STIFFER SPAR TEST POINTS 


Test points are shown for the stiffer spar, bare wing configuration. 
The testing procedure is described below. 

1. After tunnel was warmed and wind off zero readings were taken, the 
tunnel speed and density were increased to the desired values. 

2. The model was positioned at the desired angle of attack and steady 
state data was obtained. 

3. The wing was oscillated in pitch at 2, 4, 8, and 16 Hertz and un- 
steady measurements were obtained. 

4. The model was positioned at two more steady state angles of attack 
and oscillated at 2, 4, 8, and 16 Hertz. Steady and unsteady data 
were measured for each of these conditions. 

5. Tunnel speed was increased for testing at other Mach numbers for 
the same tunnel density. 

6. Upon reaching Mach 0.95 or flutter, the tunnel speed was decreased 
and Freon was pumped in to increase tunnel density to the next 
desired value. 

7. Testing resumed along another constant density line. 



0 0.5 1.0 

MACH NUMBER - M 
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NOMINAL SPAR BARE WING TEST POINTS 


Test points are shown for the nominal stiffness spar, bare wing config- 
uration. The flutter boundary for this configuration is also shown. The 
test procedure was identical to that for the stiff spar, bare wing configura- 
tion except forced oscillatory data were measured at fewer than three angles 
of attack for most test points. 
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NOMINAL SPAR WITH ENGINES TEST POINTS 


Test points are shown for the nominal stiffness spar, with engines 
configuration. The flutter boundary for this configuration is also shown. 
The test procedure was identical to that for the nominal stiffness spar, 
bare wing configuration. 
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STEADY STATE DELTA PRESSURES - INBOARD SECTION 


This is a composite plot showing how the chordwise delta pressure dis- 
tribution changes with Mach number and tunnel density. 
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STEADY STATE DELTA PRESSURES - OUTBOARD SECTION 


This is a composite plot showing how the chordwise delta pressure dis- 
tribution changes with Mach number and tunnel density. The differences 
in chordwise delta pressure distributions between the inboard and outboard 
sections is due to the difference in local angle of attack caused by the 
jig twist and the flexibility of the wing. 
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TIME HISTORY PLOTS 


This figure shows time history traces for a sample of data channels 
for both a forced oscillation case and for oscillations during flutter. 
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ROOT TORSION STRAIN GAGE - TIME HISTORY DURING FLUTTER 


This figure shows a time history plot of the root torsion strain gage 
during a flutter case. During the first 2^ seconds of this plot the model's 
deflection is increasing from flutter onset. After 2j seconds, the tunnel 
velocity was decreased by about 10% to keep the model from breaking up. 
The rest of the plot shows the model's response gradually decreasing at 
the lower tunnel speed. 


WEAKER SPAR WITH DUMMY ENGINES 

MACH = 0.67 FREQUENCY = 10 Hz 

DENSITY = 0.00172 SLUGS/FT 3 ROOT ANGLE OF ATTACK = 1.0° 
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TIME HISTORIES - BEATING FLUTTER 


Aerodynamic data was also measured while the model was beating in and 
out of flutter. This figure shows a sample of data channels during this 
beating phenomenon. 
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TIME HISTORY CHORDWISE DELTA PRESSURE OSCILLATIONS DURING FLUTTER 


Moving from the top of this figure to the bottom, one cycle of forced 
response data is shown for both the inboard and the outboard sections. 
Nine instantaneous "snapshots" are shown to depict how the chordwise delta 
pressure and airfoil position change with time. In each "snapshot" the 
top line is a bargraph of the delta pressure measurements (the leading edge 
is to the left and the trailing edge is to the right). The lower line in 
each "snapshot" depicts the unsteady airfoil position. 


FREQUENCY = 10 Hz DENSITY = 0.00172 SLUGS/FT. 3 
ROOT ANGLE OF ATTACK - +1.0° MACH = 0.67 


INBOARD OUTBOARD 

1/4 IN. - 0.5 psi 1/4 IN. = 0.5 psi 
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REPRESENTATIVE FORCED RESPONSE AERODYNAMIC COEFFICIENTS 

From the instantaneous chordwise delta pressure distribution, the lift 
and the lift and moment coefficients were calculated. These "instantaneous" 
coefficients were plotted in this figure versus the section’s "instantaneous" 
angle of attack. For the forced response case, the hysteresis moves in 
a counter-clockwise rotation indicating that energy is being put into the 
airstream by the airfoil. 





REPRESENTATIVE AERODYNAMIC COEFFICIENTS DURING FLUTTER 

This is the same type of data from oscillations during flutter. Note 
that the hysteresis is moving in a clockwise rotation indicating that energy 
is being extracted from the airstream by the airfoil. 
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REPRESENTATIVE MULTICYCLE AERODYNAMIC DATA 

Three cycles of data have been plotted to show the repeatability of 
the data. The case shown is during forced oscillations. During a divergent 
flutter case, the model's amplitude is building so the plotted data would 
also increase in magnitude. 




SECTION ANCLE OF ATTACK - (DECREES) 


563 











CHORDWISE DELTA PRESSURE COEFFICIENTS 


Unsteady chordwise delta pressures are presented in real/imaginary 

plots for representative measurements during flutter. 
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UNSTEADY LIFT COEFFICIENT - 49.6% SPAN 

Cjf o? * s are plotted for steady state, 2, 4, 8, & 16 Hertz at Ofroot=~l*5 
degrees. Magnitude/phase plots are shown for measurements at 49.6% span. 
Plots are made versus reduced frequency -k. 
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UNSTEADY LIFT COEFFICIENTS -82.1% SPAN 


Cjf^'s are plotted for steady state, 2, 4, 8, & 16 Hertz at o , roo t=-1.5 
degrees. Magnitude/phase plots are shown for measurements at 82.1% span. 
Plots are made versus reduced frequency -k. 
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UNSTEADY LIFT COEFFICIENTS 
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Cj^'s are plotted for steady state, 2, 4, 8, & 16 Hertz at aroot= _ l-5 
degrees. Real/imaginary plots are shown for measurements at 49.6% and 82.1% 
spans. Plots are made versus reduced frequency -k. 
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SUMMARY 


1. Present flutter analysis methods do not accurately predict the 
flutter speeds in the transonic flow region for wings with super- 
critical airfoils. 

2. Aerodynamic programs using CFD methods are being developed, but 
these programs need to be verified before they can be used with 
confidence. 

3. A wind tunnel test was performed to obtain all types of data neces- 
sary for correlating with CFD programs to validate them for use 
on high aspect ratio wings. The data include steady state and 
unsteady aerodynamic measurements on a nominal stiffness wing and 
a wing four times that stiffness. There is data during forced oscil- 
lations and during flutter at several angles-of-attack, Mach numbers, 
and tunnel densities. 

4. The test data is being compiled and will be published in a NASA 
report. Data will also be available through NASA on magnetic tape. 

5. The data is intended to be used for correlating with and verifying 
CFD aerodynamic programs. 


• IMPROVED TRANSONIC FLUTTER ANALYSES NEEDED 

• COMPUTATIONAL FLUID DYNAMIC CODES 

• PRESSURE / FLUTTER MODEL TEST CONDUCTED 

• COMPILATION OF TEST DATA 

• CORRELATION OF TEST DATA WITH CFD CODES 
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Aeroelastic Stability and Response of ORIGINAL PAGE IS 
Propulsion Systems (ASTROP) 0 p pQQR QUALITY 


The turbomachinery aeroelastic effort at NASA Lewis Research Center is 
focused on unstalled and stalled flutter, forced response, and whirl flutter 
of both single-rotation and counter-rotation propfans (figs. 1 & 2). It also 
includes forced response of the space shuttle main engine (SSME) turbopump 
blades (fig. 1). Because of certain unique features of propfans and the SSME 
turbopump blades, it is not possible to directly use the existing aeroelastic 
technology of conventional propellers, turbofans or helicopters. Therefore, 
reliable aeroelastic stability and response analysis methods for these 
propulsion systems must be developed. The development of these methods for 
propfans requires specific basic technology disciplines, such as 2D and 3D 
steady and unsteady (unstalled and stalled) aerodynamic theories in subsonic, 
transonic and supersonic flow regimes; modeling of composite blades; geometric 
nonlinear effects; and passive or, active control of flutter and response. 

These methods for propfans are incorporated in a computer program ASTROP 
(fig. 3). The program has flexibility such that new and future models in 
basic disciplines can be easily implemented. The forced response analysis 
method for turbine blades will be discussed later. 


UNSTEADY AERODYNAMICS AND AEROELASTICITY 


SR LARGE SCALE 
AERODYNAMICS AND FLUTTER 
SR 7 L 



SR-3 RESEARCH MODEL 
FLUTTER, FORCED 
RESPONSE, AND MISTUNING 



CR PROPFAN WIND TUNNEL MODEL 
AERODYNAMICS AND FLUTTER 



SSME TURBOPUMP BLADE 
FORCED RESPONSE 



Figure 1 
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UNSTEADY AERODYNAMIC AND AEROELASTICITY 
FLIGHT TEST PROGRAM 


CR PROPFAN 


SR PROPFAN 
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Figure 2 


ASTROP 



Figure 3 
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Application of the ASTROP Code to Investigate Flutter 
of a Composite SR Propfan Model 

One version of the code called ASTR0P3 (ref. 1) uses three-dimensional 
subsonic steady and unsteady cascade aerodynamics (ref. 2) and NASTRAN (ref. 3) 
finite element model to represent the blade structure. The equivalent aniso- 
tropic material properties for each finite element are generated by using a 
preprocessor code COBSTRAN (ref. 4). The effect of centrifugal loads and 
steady-state airloads on the steady-state geometry of a composite wind tunnel 
model (SR3C-X2) blade is shown in figure 4(a). The aerodynamic cascade 
effects (or the effect of number of blades) on the eigenvalues are shown in 
figure 4(b). Both centrifugal loads and aerodynamic loads untwist the blades 
and this untwist increases with rotational speed. It is evident from fig- 
ure 4(b) that the number of blades or the cascade effect is very significant 
on the real part of the eigenvalue and hence on stability. 

EFFECT OF CENTRIFUGAL AND AERODYNAMIC LOADS ON 
BLADE PITCH ANGLE 



(a) 


EFFECT OF NUMBER OF BLADES ON EIGENVALUES 



Figure 4 (b) 


574 



Comparison of Measured and Calculated Flutter 
Boundaries for the SR3C-X2 Propfan Model 

Theoretical flutter results obtained from the ASTR0P3 code have been 
correlated in figures 5(a) and (b) with flutter data of a wind tunnel propfan 
model (ref. 5), SR3C-X2 , with composite blades. Theoretical results include 
the effects of centrifugal loads and steady-state airloads. The theory does 
reasonably well in predicting flutter speeds and slopes of the boundaries. 
However, the difference between the calculated and measured flutter Mach num- 
bers for the four-blade case is greater than for the eight-blade case. This 
implies that the theory may be overcorrecting for aerodynamic cascade effects 
for four blades. Calculated interblade phase '-angles at flutter (not shown) 
also compared well with measured values. However, calculated flutter fre- 
quencies were about 8% higher than measured. 


4 BLADES 


8 BLADES 



(a) 


Figure 5 


(b) 
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Evaluation of Two-Dimensional Unsteady Aero 
for Propfan Flutter Prediction 

Actually the ASTROP code was started with two-dimensional unsteady 
aerodynamic theory (ref. 6) by correcting for blade sweep (ref. 1). The 
version of the code which uses blade normal modes and two-dimensional unsteady 
aero theory in a stripwise manner is ASTR0P2 (figure 3). The ASTR0P3 version 
uses three dimensional unsteady aero theory. To assess the validity of 
two-dimensional aerodynamic theory and the associated sweep correction, the 
real part of the eigenvalue of the critical mode calculated by using both 
ASTR0P2 and ASTR0P3 are compared in figure 6. Also included in this figure is 
the measured flutter mach number. Evidently, the two-dimensional theory is 
less accurate than three-dimensional theory in predicting flutter Mach number 
for this case. Correlative studies (not shown) of measured and calculated 
flutter boundaries were also conducted by varying Mach number, blade sweep, 
rotational speed, and blade setting angle. The correlation varied from poor 
to good. In some cases the expected conservative nature of the 
two-dimensional theory did not prevail, possibly because of the arbitrary 
nature of the reference line which is employed in the strip-method, and the 
associated sweep correction. 

COMPARISON OF 2-D AND 3-D UNSTEADY AERO FOR 
PROPFAN FLUTTER PREDICTION 
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Propfan Blade Mistuning Models 

Understanding the effects of blade mistuning on vibration, flutter and 
forced response of turbomachinery rotors is a current research topic because 
mistuning affects flutter and response behavior. An analytical and 
experimental investigation of mistuning in propfan flutter was conducted in 
ref. 7. A schematic for an eight-bladed mistuned rotor in formulating the 
analytical model and blade ply directions used in constructing the wind tunnel 
model are shown in figures 7 (a) and (b), respectively. The analytical model, 
which is more general than the wind tunnel model, is based on normal modes of a 
rotating composite blade and subsonic unsteady lifting surface aerodynamic 
theory. The natural frequencies and mode shapes of the SR3C-X2 and -3 model 
blades differ because of the ply angle variations between the blades. The 
first mode frequencies of both the blades are very close and were insensitive 
to ply angles. However, the average second mode frequency of the SR3C-3 blade 
is about 12 percent higher than that of the -X2 blade. More details can be 
found in ref. 7. 

BLADE GROUP SCHEMATIC FOR AN EIGHT-BLADED ROTOR BLADE PLY DIRECTIONS 



NG = i) 
NB - 2 


(a) 


Figure 7 
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Comparison of Eigenvalues of Tuned 
and Mistuned Propfan Models 

To illustrate the effect of mistuning (which is partly aerodynamic and 
structural) the calculated real and imaginary parts of eigenvalues of the 
SR3C-X2 (8-bladed tuned rotor), SR3C-3 (8-bladed tuned rotor), and mixed 
(mistuned) rotor were compared in figure 8. The mistuning is due to the dif- 
ferences in blade steady-state geometry, frequencies, and mode shapes. The 
eigenvalues are for all interblade phase angles of the mode with lowest 
damping. The calculations were performed by treating the SR3C-X2 and -3 
rotors as tuned and the mixed rotor as an idealized alternately mistuned 
rotor--four identical blade pairs with two different blades in each pair. 
Comparison of root loci indicates that the area of the approximate ellipse 
for the SR3C-X2 is greater than that of SR3C-3, indicating a stronger aero- 
dynamic coupling between the blades of the SR3C-X2 rotor. The difference in 
stability of the tuned rotors is due to the difference in blade stiffness and 
mode shapes because of the different ply angles of the blades. The results 
also show that mixing the blades significantly affected the eigenvalues and 
resulted in a rotor with a greater damping than the lowest damped mode of 
either tuned rotor. 


ROTOR 
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REAL PART OF EIGENVALUE 
Figure 8 


Comparison of Measured and Calculated Flutter Boundaries 
for Tuned and Mistuned Propfan Models 

Measured and calculated results for the tuned rotor SR3C-X2 and mistuned 
rotor SR3C-X2/SR3C-3 are compared in figure 9. The calculations for each 
rotor were made with the calculated modes and frequencies, except that the 
measured second mode frequency was substituted for the calculated one. The 
calculated flutter Mach numbers for the SR3C-X2 are less than the measured 
ones for all rotational speeds. The agreement would be better if the effects 
of steady airloads and structural damping were included in calculations (see 
ref. 1 for detailed discussion). The agreement of the mixed rotor is better, 
but would become unconservative if steady airloads and structural damping were 
included in the theory. However, the overall agreement between theory and 
experiment is more than satisfactory. For additional details and results see 
ref. 7. 



Figure 9 


Comparison of Measured and Calculated Vibratory 
Stress Amplitudes of a Propfan Model 

A new feature of the ASTR0P3 code under development is the capability to 
perform a modal forced response vibration analysis of aerodynamical ly excited 
propfans. Figure 10 depicts a single-rotation, advanced propfan wind tunnel 
model (SR5, 10 metallic blades), ref. 8, operating in a generally uniform, 
steady inflow field, inclined at a small angle with respect to the axis of 
rotation. Although the absolute inflow field is constant, rotation of the 
propfan results in velocities with oscillatory components relative to the 
blades. Under such conditions, ASTR0P3 is able to determine the oscillatory 
loading distributions over the propfan blades at various excitation fre- 
quencies and calculate the vibratory displacements and stresses of the prop- 
fan. The table shows measured and preliminary calculated one per rev 
vibratory stress amplitudes for the SR5 blade. Also included in the table 
are the calculated results from ref . 9 by using a 2-D unsteady aerodynamic 
theory. Comparison shows that 2'-D results are better than 3-D results. The 
reason for this difference is being investigated. 
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Stall Flutter Analysis Methods 


The third feature of the ASTROP code is a stall flutter analysis which is 
in ASTROPS. Under take-off conditions, the propfan blades operate at high 
angles of attack and have the potential to stall flutter, triggered by 
separated flow during part of every cycle of oscillation. Stall flutter 
speeds are very low and the forces due to vibration at the stall condition 
(dynamic stall) are an order of magnitude high compared to forces in separated 
flow with no vibration. Prediction of forces during dynamic stall has been a 
continuing research effort. Some prediction methods are reviewed in ref. 10, 
and their classification is shown in figure 11. The Navier-Stokes solvers 
(N.S.S.), vortex methods, and the zonal methods attempt to solve the fluid 
mechanics equations in their fundamental form by numerical techniques with 
varying degrees of simplifications and assumptions. These models require a 
significant amount of computer time and therefore are not suitable for routine 
aeroelastic analysis. In semi-empirical models an analytical approximation is 
attempted to approximately reproduce measurements for example, by way of 
analytical curve fit to wind tunnel data. The semi-empirical models take less 
computer time to solve and can be used in a routine aeroelastic analysis. 



Figure 11 
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Comparison of Measured and Calculated Stall Flutter Results 
of an Unswept Propfan Model 

Recently two semi-empirical dynamic stall models, designated as model A 
and model B, were applied to an unswept propfan model (SR-2, 8 metallic 
blades). The calculated logarithmic decrement of third mode response as a 
function of blade pitch angle is shown in figure 12. Also included in the 
figure is the blade pitch angle at which the blade is unstable in a wind 
tunnel experiment. Model A (ref. 11) incorporates the unsteady effects in 
stall using only one stall parameter that relates the dynamic stall angle and 
the non-dimensional rate of angle-of-attack. The stall parameter is given as 
a function of Mach number and airfoil thickness to chord ratio. Model B (ref. 
12) is a synthesized data method to dynamic stall modeling. An analytical 
curve is fitted for the wind tunnel data obtained from oscillating airfoil 
tests. The empirical parameters in the model are obtained from this fit. 
However, experimental dynamic data is not available for propfan airfoil 
sections (16 series). Therefore, in implementing model B for propfan 
application, the dynamic data corresponding to NLR-1 airfoil was used, even 
though the airfoil geometries are different. However, for the case studied 
here, the Mach number range for which the data available for the NLR-1 airfoil 
corresponds to the helical Mach number (at zero freestream Mach number) due to 
rotation of the propfan model. In spite of the differences in dynamic data, 
the empirical models chosen predicted a qualitative stall flutter behavior for 
the case studied. Both the models predicted that the stall flutter response 
occurred in third mode as was found in the experiment at a rotor speed of 8500 
rpm. However, the calculated blade pitch angles at which the stall flutter 
occurred is lower than that of the experiment (30°), model B predicted a 
closer value (28.25°) compared to that predicted by model A (25°). The 
calculated frequency at stall flutter condition is about 10% higher than the 
experimental data. 



Figure 12 
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Stall Flutter Analysis with Navier-Stokes Code 

A solution procedure is described for determining the two-dimensional, 
one- or two-degree-of-freedom flutter characteristics of arbitrary airfoils at 
large angles-of-attack (see ref. 12). The same procedure is used to predict 
stall flutter including separated flow. This procedure requires a simul- 
taneous integration in time of the solid and fluid equations of motion. The 
fluid equations are the unsteady compressible Navier-Stokes equations, solved 
in a body-fitted moving coordinate system using an approximate factorization 
scheme. The solid equations are integrated in time using an Euler implicit 
scheme. Several special cases, figures 13-15, are presented to demonstrate 
the capability of this scheme to predict transonic flutter and stall flutter 
with large separated flow. 

The first case is shown to illustrate the capability of the present solver 
to predict the highly separated flows. The aerodynamic coefficients of an 
NACA 0012 airfoil oscillating in pitch at large angle-of-attack is shown and 
compared with experiment in figure 13. The mean angle and amplitude of oscil- 
lation was 15 degrees and 10 degrees respectively. The reduced frequency 
based on semi-chord was 0.151. The freestream Mach number and Reynolds num- 
ber were 0.283 and 3.45 million respectively. It is seen from figure 13 that 
the Navier-Stokes solver produces lift, drag, and moment coefficients which 
are in a reasonable agreement with the measured ones. The fact that the flow 
solver is able to capture much of the dynamic stall flow features increases 
the confidence in the capability of this code for stall flutter predictions. 
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STALL FLUTTER ANALYSIS WITH NAVIER-STOKES 
CODE: a = 15 - 10 cos (fit) 
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Comparison of Flutter Speeds from Euler, UTRANS2, 
and LTRAN2 Codes 

A second special case considered for validation of the Navier-Stokes 
solver is its Euler version to calculate transonic flutter speed. Predicted 
transonic flutter speeds at various mass-to-air ratios for NACA 64006 airfoil 
oscillating in pitch and plunge at Mach number 0.85 are shown in fiqure 14. 
Several other validation cases are reported in ref. 13. The results from 
UTRANS2 (ref. 14) and LTRAN2 (ref. 15) are also included. Very good agreement 
between present Euler and UTRANS2 code results is found. A qualitative 
agreement between present and LTRAN2 results is found, too. 



Figure 14 
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Plunging and Pitching Stall Response 


The third case considered with the Navier-Stokes code is the prediction 
of flutter at large mean angle-of-attack, including flow separations. The 
time response of plunging and pitching displacements and lift and moment 
coefficients of an NACA 0012 airfoil is shown in figure 15. The airfoil was 
initially subjected to a sinusoidal pitching oscillation from 5 to 25 
degrees. During the downstroke, around 23.8 degrees, the airfoil was released 
and was allowed to follow a pitching and plunging motion. The dimensionless 
speed is varied from 4 to 8. The response of the airfoil is stable when the 
speed is 4 and is unstable when the speed increases to 8. It was found that 
the growing response is induced by the separated flow over the airfoil at 
large angle-of-attack. 
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Full Potential Unsteady (Cascade) Aero Model 

A compressible, unsteady, full potential, finite difference code is being 
developed for modeling 2D/3D flow through single rotation propfans and other 
turbomachinery rotors. The procedure introduces a deforming grid with a uni- 
form shear mesh. The numerical scheme is based on finite volume and implicit 
time marching technique. The 2-D code is vectorized and verified by applying 
it to several special cases. Two such cases are shown in figure 16. For 
comparison, the results from refs. 16 and 17 are also included. Even though a 
very coarse grid is used in the present calculations, the agreement between 
the present results and those of refs. 16 and 17 is very good. Validation of 
3-D code is in progress. 
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Figure 16 
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CR Propfan 3-D Steady and Unsteady Aero Model 

An incompressible, steady and unsteady aerodynamic model is being devel- 
oped in ref. 18 for modeling flow through counter rotation propfans. The model 
■is based on time domain solution in conjunction with panel method. This model 
is applicable for calculating performance and stability of both single and 
counter rotation propfans including interaction from wing. The code is being 
verified by applying to several special cases. One such case is shown in fig- 
ure 17 in which present results are compared with the corresponding ones in 
ref. 19. See ref. 18 for additional validations and for details. This code 
will be extended to compressible flow, and, then, will be merged with ASTROP 
code structural modules to predict flutter of counter rotating propfans. 
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Propfan Wind Tunnel Flutter Models 

Figure 18 shows propfan blade wind tunnel models that have had flutter. 
These blades are not aeroelastical ly scaled models and were made for 
aerodynamic performance tests. However, the SR3C-X2 and -3 models were 
specifically designed for flutter and forced response experiments, 
respectively. One single rotation model had stall flutter, SR-2. The other 
two had unstalled flutter SR3C-X2 and SR-5. The flutter data from these 
models has been used to verify the analysis methods discussed earlier. The 
three counter-rotation models shown have had unstalled flutter. The 
correlation of this data with analysis is in progress. 
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Forced Response of SSME Turbopump Blades 


The state-of-the-art in forced response analysis of turbomachinery 
blading is to simply calculate the blade natural frequencies and try to avoid 
known forcing frequencies (Campbell Diagram).' Forced response calculations 
are not attempted. This can lead to unexpected blade cracking. The objective 
of this research is to develop a forced response prediction method for 
turbopump blades. The flow chart for this method is shown in figure 19. The 
development will proceed in three parallel, integrated tasks. The first task 
continues existing in-house research to develop a model (M-Stage) of the 3D, 
time-averaged, flow field within a passage of a blade row embedded in a 
multi-stage machine. This model identifies the distorted (i.e. 
non-axisymmetric) flow field generated by neighboring blade rows. This 
information serves as input to Task 2. The second task will develop a model 
(LINPOT) to predict the unsteady aerodynamic loads generated by the flow 
distortion. This model will consist of an unsteady, linearized, potential 
flow solver, and a linearized, convected gust solver. The model will be 
applicable to thick, highly cambered turbine blades. The third task continues 
existing in-house research to develop a model (FREPS) for integrated forced 
response predictions. This' model will integrate the M-STAGE model of task 1 
with the LINPOT model of task 2 and a structural dynamic model. Two 
structural dynamic models will be used. Initially, a simplified two 
degree-of-f reedom blade model will be incorporated. This will be followed by 
a complete modal blade model. The result of this research will be a system to 
calculate the forced response of a turbopump blade embedded in a multi-stage 
turbine. The benefit will be a marked reduction in occurrences of unexpected 
blade cracking. This system will also be applicable to blading in 
aeronautical propulsion systems. 



Figure 19 
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20 Unsteady, Viscous Cascade Aero Model 

A compressible, unsteady, full Navier-Stokes, finite difference code has 
been developed for modeling transonic flow through two-dimensional, 
oscillating cascades. The procedure introduces a deforming grid (fig. 20) 
technique to capture the motion of the airfoils. The use of a deforming grid 
is convenient for treatment of the outer boundary conditions since the outer 
boundary can be fixed in space, while the inner boundary moves with the blade 
motion. The code is an extension of the isolated airfoil code developed in 
reference 20. More results validating the deforming grid technique are 
presented in reference 21. 

2D UNSTEADY VISCOUS CASCADE AERO MODEL 



2D UNSTEADY VISCOUS CASCADE AERO MODEL 



Figure 20 
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Pressure Coefficients for an NACA 0012 Cascade 
with Viscous Effects 

The deforming grid technique has been used to predict the load histories 
for a NACA 0012 cascade with zero inter-blade phase angle and zero stagger. 

Two flow conditions were selected to investigate both subsonic and transonic 
flow. The cascade has a gap to chord ratio of one, M = 0.60 and 0.67, Re = 
3.21 million, ct^ = 0.0 degrees, pitching t2.0 degrees, and k = 0.20 
(reduced frequency based on semi-chord). A Fourier transform on the pressure 
coefficient distribution was done for the first harmonics. The results are 
shown in figure 21. Future work will investigate non-zero inter-blade phase 
angles and will compare predictions with experimental data from the NASA Lewis 
Transonic Oscillating Cascade Facility. 
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3-D Unsteady Euler Analysis 

The three dimensional unsteady Euler analysis for an isolated wing 
developed in ref. 22 has been recently extended to propfans. This extended 
version of the code is being utilized to study the effect that propfan 
angle-of-attack has on the unsteady blade loading of a single rotation propfan 
design. The code is capable of modeling the complete propfan configuration. 
The program will be used to predict the unsteady loading on the propfan 
recently tested in a two-bladed configuration as part of the Large-Scale 
Advanced Propfan (LAP) program. Part of the objective of this test program 
was to obtain detailed steady and unsteady blade surface pressure measurements 
for benchmarking computer models. Presently this computer program assumes 
that the blades are rigid. It is planned to look at the formulation and 
coding necessary to allow the blades to respond to the unsteady loading thus 
allowing the program to be used in aeroelastic forced response predictions. 
Sample pressure contours on blades of a propfan are shown in figure 22. 



Figure 22 
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2-D Unsteady Perturbation Analysis for Cascades 


In order to study the effect that blade sweep has on the flutter behavior 
of a cascade of airfoils operating in the transonic flow regime, the 
linearized unsteady analysis developed in ref. 17 is being utilized. This 
analysis predicts the unsteady loading resulting from small amplitude harmonic 
motion of the blades in a two-dimensional cascade operating in an inviscid 
subsonic or transonic flow. The unsteady potential is obtained from a 
perturbation analysis applied to the steady flow solution. Thus, the unsteady 
analysis is able to include the effects of finite mean loading on the unsteady 
response. At LeRC, the transonic potential code developed in reference 23 is 
utilized in calculating the steady flow field. Sample unsteady pressures 
calculated for a cascade of NACA 0012 airfoils (at Mach number 0.6, stagger 
angle 45° and mean incidence angle 9°) by using the combined code are shown in 
figure 23. The combination of these steady and unsteady programs allows for 
the prediction of the flutter behavior of fan and propfan designs which 
include the effect of realistic reduced frequencies and blade geometries. The 
resulting computer program will be benchmarked against experimental cascade 
data and then applied to study the effect that blade sweep has on the flutter 
behavior of a cascade of airfoils. 
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Unsteady Swept Cascade Experiments 

The effect of blade sweep on the flutter behavior of a cascade of 
airfoils is being studied in the transonic oscillating cascade facility, 
figure 24. This study is being conducted to determine if classical sweep 
corrections applicable for fixed wing flows are valid for oscillating airfoils 
in a torsional motion while maintaining a selected interblade phase angle 
between adjacent blades. The unsteady loading is determined by a number of 
blade mounted high response pressure transducers. The initial phase of 
testing will involve the use of unswept airfoils in order to provide a 
baseline set of data for benchmarking the computer programs to be used in this 
study. The swept airfoils will then be installed and a series of tests will 
be run to determine the effect of incidence angle, Mach number, reduced 
frequency and interblade phase angle on the flutter behavior of the swept 
cascade. 



Figure 24 
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Three-Dimensional Gust Model for a Propeller Blade 

A propeller blade (figure 25) rotating in a nonuniform upstream flow 
encounters an unsteady flow field, even when the nonuniform upstream flow is 
steady. For straight bladed propellers, the unsteady flow of the propeller 
can be approximated by a two-dimensional wing in a three-dimensional gust as 
shown in figure 25. For small amplitude disturbances, the unsteady flow field 
may be obtained as a perturbation about the underlying steady flow. The 
governing equation is a linear, nonconstant coefficient, inhomogenous, 
convective wave equation, see refs. 24 and 25, A finite difference scheme is 
used to solve for the perturbation potential. Some sample results are shown 
in figure 26. 



Figure 25 
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Variation of Unsteady Lift Coefficient 
of 1 2 % Thick, Symmetric Joukowski 
Airfoil in a Transverse Gust 

The governing wave equation for the model described in figure 25 is 
solved for perturbation potential by using a finite differencing scheme. For 
a symmetric airfoil in a transverse gust the real and imaginary parts of the 
lift at Mach number 0.6 and with reduced frequency as parameters are shown in 
figure 26. Also included in the figure is the corresponding curve for the 
flat plate. Comparing the flat plate and 12% thick airfoil results, it is 
observed that the thickness effects on the lift are more significant at low 
reduced frequencies. Similar comparisons (not shown) are also made at 
different Mach numbers, and it was found that the thickness effects on lift are 
more significant at higher Mach numbers. 
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Numerical Simulation of Flow Through 
Counter Rotating Propfans - Average 
Passage Flow Model 

A new analytical model, referred to as "average passage flow model," is 
being developed in reference 26 for simulating flows through counter rotating 
propfans. It describes the three-dimensional time-average flow field within a 
typical passage of a blade run in a multiblade run configuration. The model 
has been used to examine the flow field generated by a counter rotating 
propfan configuration (UDF). For example, the pressure field radiated by the 
aft fan is shown in figure 27. The pressure field is color-coded with a 
spectrum ranging from blue (low pressure) to green to yellow to red (high 
pressure). Tne boundary between green (lower pressure) and yellow-orange 
(higher pressure) regions is the footprint of the aft propfan trai ling-edge 
shock. The base of the shock lies at approximately three-quarters of the 
span. From this point it appears to spiral outward beyond the tip of each 
bl ade. 



Figure 27 
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Application of "Average Passage Flow Model" 
for CR Propfan Noise Prediction 

The average passage flow model developed in ref. 26 was merged with an 
aeroacoustic prediction model for CR propfans developed by Dr. F. Farassat of 
NASA Langley Research Center. This merger permits the simultaneous evaluation 
of aerodynamic performance and radiated sound levels. Figure 28 shows a 
comparison of the predicted sound levels with corresponding measured data of 
the CR scaled model by Dr. J. H. Dittman of NASA Lewis Research Center. The 
correlation between theory and experiment is excellent. 
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Figure 28 
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LeRC Groups Involved in Unsteady Aerodynamics and 
Aeroelasticity 

The groups in the unsteady aerodynamics and aeroelasticity effort at LeRC 
are shown in figure 29. Also listed in the figure are the names of NASA LeRC 
employees. Support Service Contractors and Grantees, who contributed to the 
research effort described in the paper. 
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Figure 29 
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LeRC Unsteady Aerodynamics and Aeroelasticity Effort 

The elements of the overall research effort in the subject area are shown 
in figure 30. These elements cover the development of unsteady aerodynamic 
models, aeroelastic models (for flutter, forced response and optimization), 
associated computer programs, and wind tunnel flutter experiments. 
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INTRODUCTION 


The requirements placed on an unsteady aerodynamic theory intended for turbomachinery aero- 
elastic or aeroacoustic applications will be discussed along with a brief description of the various 
theoretical models that are available to address these requirements. The major emphasis is placed on 
the description of a linearized inviscid theory which fully accounts for the affects of a nonuniform 
mean or steady flow on unsteady aerodynamic response. Although this linearization has been 
developed primarily for blade flutter prediction, more general equations will be presented which 
account for unsteady excitations due to incident external aerodynamic disturbances as well as those 
due to prescribed blade motions: In this presentation we will focus on the motivation for this 

linearized unsteady aerodynamic theory, outline its physical and mathematical formulation and 
present examples to illustrate the status of numerical solution procedures and several effects of 
mean-flow nonuniformity on unsteady aerodynamic response. This presentation is based on a paper 
of the same title which is published in full in the Proceedings of the Tenth U. S. National Congress 
of Applied Mechanics (ref 1). 


• Linearized unsteady aerodynamic analysis 

• Real blade geometry 

• Mean blade loading 

• Shock phenomena 

• Design applications 

• Aeroelastic 

e Blade flutter and forced vibration 

• Aeroacoustic 

• Noise generation, transmission 
and reflection 
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MAJOR ASSUMPTIONS 


The development of theoretical models to predict unsteady flows through turbomachines is a 
formidable task. The analyst is confronted with determining the time-dependent, three-dimensional 
flow of a viscous compressible fluid through a geometric configuration of enormous complexity. This 
task has required the introduction of a considerable number of simplifying assumptions to make 
the problem mathematically tractable and to render the resulting solutions useful to designers. 
For the most part, the theoretical formulations that have been developed to predict the unsteady 
aerodynamic phenomena associated with blade flutter or forced vibration consider the blades of 
an isolated two-dimensional cascade, neglect viscous effects at the outset and regard unsteady 
fluctuations to be of sufficiently small amplitude so that a linearized treatment of the unsteady 
perturbation is justified. In addition the resulting two-dimensional inviscid flow is assumed to 
remain attached to the blade surfaces, the mean flow is assumed to be at most a small irrotational 
steady perturbation from a uniform stream at the cascade inlet, and any shocks that might occur 
are assumed to be of weak to moderate strength and have small curvature. 


• Isolated blade row 

• Two-dimensional inviscid (Re-^oo ) flow 

• Small-amplitude periodic unsteady 
excitations 

• Attached flow 

• Irrotational mean flow at Inlet: 

7= v m00 + v? 

• Weak shocks 
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REQUIREMENTS 


In general, the unsteady aerodynamic analyses intended for turbomachinery aeroelastic appli- 
cations must be applicable to fan, compressor and turbine cascades, to subsonic, transonic and 
supersonic Mach numbers and to moderate through high frequency structural and external aero- 
dynamic excitations. Then, to determine the aeroelastic and aeroacoustic characteristics of the 
blading such analyses must be capable of predicting the unsteady loads acting on the blades and 
the amplitude and wave numbers of the acoustic waves which carry energy away from the blade 
row and the entropic and vortical fluctuations which are convected downstream. These responses 
arise from the various sources of unsteady excitation including prescribed blade motions, variations 
in total temperature and pressure (“entropy and vorticity waves”) at inlet, and variations in static 
pressure (acoustic waves) at inlet and exit. For blade flutter applications it is only necessary to 
predict the unsteady loads acting on the blades as a result of prescribed blade motions; for forced 
response applications the unsteady blade loads due to incident entropic, vortical and acoustic dis- 
turbances are also required. Finally, for aeroacoustic applications the parameters associated with 
far-field acoustic responses must be determined. 


• Fan, compressor and turbine cascades 

• Subsonic, transonic, supersonic Mach numbers 

• Moderate to high excitation frequencies 

• Response predictions 

• On blades: surface pressures, global unsteady airloads 

• Far field: outward propagating acoustic waves 
vorticity and entropy variations downstream 

• Prescribed excitations 

• Blade motions (flutter) 

• External aerodynamic disturbances (forced vibration) 

• Vortical and entropic disturbances at inlet 

• Acoustic disturbances at inlet and exit 


i 
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TRANSONIC CASCADE: < M_ 0 0 < 1 


A representative cascade configuration is shown in the figure below which depicts a two-dimen- 
sional section of a transonic compressor rotor ( M < M_ 00 < 1. The cascade stagger angle is 
denoted by 0 and the blade spacing by G. In the absence of unsteady excitation the blades are 
identical in shape, equally spaced and their chord lines are oriented at the same angle, 0, relative 
to the axial flow direction. The inlet and exit free-stream flows are described by the velocity vectors 
V^poo. The free-stream flow angles measured relative to the axial-flow (or £-) direction are denoted 
by For the configuration illustrated the inlet and exit conditions are such that normal shocks 

(Sh) emanate from the blade suction surfaces and vortex wakes (W) emanate from the blade trailing 
edges and extend downstream. 
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UNSTEADY EXCITATIONS 


The unsteady fluctuations in the flow arise from one or more of the following prescribed sources: 
blade motions, upstream and/or downstream acoustic disturbances which carry energy toward the 
blade row, and upstream entropic and vortical disturbances which are convected through the blade 
row by the mean flow. These excitations are assumed to be of small amplitude and periodic in time. 
The external aerodynamic excitations are also spatially periodic, while the structural excitation 
is periodic in the “circumferential” or 77-direction. For example, we consider blade motions and 
incident acoustic disturbances as described below. Here 7 ? measures the displacement of a point on 
a moving blade surface relative to its mean or steady-state position, X is a position vector, m is a 
blade number index, t is time, is the reference-blade (m = 0) complex displacement-amplitude 
vector, u> and <7 are the temporal frequency and interblade phase angle, respectively, of the unsteady 
excitation and Re{ } denotes the real part of { }. Also, p/ )TO o and ic Too are the amplitude and wave 
number, respectively, of an incident pressure fluctuation, pj(X,t), coming from far upstream (—00) 
or far downstream (+00). Note that the interblade phase angle, <7, of an incident disturbance is 
The temporal frequency and wave number of an incident vortical or entropic disturbance 
are related by = — « Too • G, but a more complicated relationship exists between u> and « ?0O for 
an incident pressure disturbance. 



Blade motions: R (X + mG, t) = Re j r*(Xj exp [i (cot + mo)] J 
Incident disturbances: p'- «> (X, t) = Re jp^ooexp [i(/c_ oo* X + cot)] | 
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TIME-DEPENDENT GOVERNING EQUATIONS 


The equations governing the fluid motion follow from the integral forms of the mass, momentum 
and energy conservation laws and the thermodynamic relations for a perfect gas. The former provide 
a coupled set of corresponding nonlinear differential equations (the Euler equations) in continuous 
regions of the flow and jump conditions at surfaces across which the inviscid flow variables are 
discontinuous, i.e., at vortex-sheet wakes and shocks. In continuous regions the energy equation can 
be replaced by the requirement that the entropy following a fluid particle must remain constant. 
In addition to the foregoing field equations and jump conditions, the attached flow assumption 
requires that the unsteady flow must be tangential to the moving blade surfaces and information 
on the uniform inlet and exit flow conditions and the incident entropic, vortical and acoustic or 
static pressure disturbances must be specified. The remaining steady and unsteady departures from 
the uniform inlet and exit conditions must be determined as part of the time-dependent solution. 
This foregoing aerodynamic problem is a formidable one as it involves a system of nonlinear time- 
dependent equations with conditions imposed on moving blade, wake and shock surfaces in which the 
instantaneous positions of the wakes and shocks must be determined as part of the solution. Because 
of these features and the prohibitive expense that would be involved in obtaining the aerodynamic 
response information needed for aeroelastic or aeroacoustic applications, the usual approach is to 
examine limiting forms of the full governing equations with the intention of providing efficient 
analyses for design applications. 


• Integral conservation laws 

• Thermodynamic relations 



• Euler equations at field points 

rsj 

• Jump conditions at moving shocks (Mf ± 0) and at 

vortex sheet boundary layers and wakes (Mf = 0) 

• Flow tangency condition at moving blade surfaces 



• Far-field behavior 
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UNSTEADY AERODYNAMIC LINEARIZATIONS 


Because of the complexity of the nonlinear time-dependent unsteady aerodynamic problem, 
linearized treatments of the unsteady flow are often considered. The major linearizations that 
have been proposed are the following: classical linearized theory, time-linearized transonic flow 
theory and the present theory in which unsteady disturbances are regarded as small relative to 
a fully nonuniform mean flow. The essential differences between these theories arise from the 
manner in which the steady flow is represented. In classical theory both steady and unsteady 
departures from a uniform stream are regarded as small and of the same order of magnitude. 
In time-linearized transonic theory steady and unsteady disturbances are regarded as small and 
very small, respectively, relative to uniform free-stream flow properties. Finally, in the present 
linearization no restriction is placed on the steady flow but the unsteady perturbations are assumed 
to be of small amplitude. The classical theory applies at the (reduced) frequencies of interest for 
turbomachinery applications, but steady flow variations have no impact on the unsteady response. 
Time-linearized transonic theory applies at Mach numbers near one and the unsteady perturbation 
depends on the steady flow, but this theory is formally restricted to low-frequency unsteady motions. 
The present theory fully includes the effects of nonuniform mean flow and applies throughout the 
Mach number and frequency range of interest for turbomachinery applications. 


• Classical theory 

P U,y,t) = P-oo + p(x,y) + Re | p(x,y)e lwt | + •• 

• Time-linearized transonic theory 

p ( *i y «t ) * P-® + P ( x,y ) + Re |p(x,y)e ia ' t i+... 

• Present theory 

P (x,y,t) = P(x,y) + Re j p(x,y)e lwt ( + . . . 
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PRESENT LINEARIZATION 


The equations governing small-amplitude unsteady departures from a nonuniform mean flow are 
determined by expanding the various flow variables in asymptotic series in e, where e is a measure of 
the amplitude of the unsteady excitation. Thus, for example, in the first equation below P(X,t) is 
the time-dependent fluid pressure, P(X) is the pressure in the steady background flow, Re{p(X)e lujt ] 
is the first-order time-dependent pressure and p(X) is its complex amplitude. In addition, Taylor 
series expansions and surface vector relations are used to refer information at a moving blade, wake 
or shock surface (<S) to the mean position of this surface (S). Equations governing the zeroth-order 
or steady and the complex amplitudes of the first-order unsteady flow properties are obtained after 
substituting the foregoing expansions into the full time-dependent governing equations, equating 
terms of like power in e and neglecting terms of higher than first order in e. It follows from the 
original assumptions that the steady background flow is governed by a full-potential boundary- 
value problem and that the complex amplitudes of the unsteady flow properties are governed by 
a system of time-independent linearized equations with variable coefficients which depend on the 
underlying mean flow. In the unsteady problem surface conditions can be imposed at the mean 
surface locations, and in both the steady and first-order unsteady problems the required solution 
domain can be limited to a single extended blade- passage region. 


• Series expansions 

P (X,t)= P(X) + Re {p(X)e ia,t }+ ... 



(P + /? .yP + ...) s 



fg* f (t s ,/? ) + 


• Full potential boundary-value problem for steady 
background flow 

• Linear variable-coefficient boundary-value problem for first- 
order unsteady flow 

• Time-independent 

• Surface conditions at mean surface locations 

• Single extended blade-passage solution domain 


I 
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THE STEADY BACKGROUND FLOW 


As a consequence of our assumptions regarding shocks and the steady flow far upstream of 
the blade row, the mean or steady background flow through the cascade will be homentropic and 
irrotational; i.e., V = where V and $ are the local steady velocity and velocity potential, 
respectively. The field equations governing the steady flow follow from the mass and momentum 
conservation laws and the isentropic relations for a perfect gas. Here, p, P , M and A and are the 
local steady density, pressure, Mach number and speed of sound propagation, respectively, and 7 is 
the specific heat ratio of the fluid. Surface conditions for the zeroth-order or steady flow apply at 
the mean positions of the blade ( B ), wake (W) and Shock (Sh) surfaces. Blade mean positions are 
prescribed, but the mean wake, i.e., the stagnation streamlines downstream of the blade row, and 
shock positions must be determined as part of the steady solution. Since, by assumption, the flow 
remains attached to the blades, a flow tangency condition applies at blade surfaces. In addition 
the steady pressure and normal velocity component must be continuous across blade wakes and 
mass and tangential momentum must be conserved at shocks. Finally, three of the far-field uniform 
velocity components, or the equivalent information, must be prescribed to completely specify the 
steady boundary-value problem. The fourth or remaining component can be determined in terms 
of the three prescribed using the integral form of the mass conservation law. 


• Field equations : V(/>VS>) s O 

p {r ~' ] ’(rMljP/2) (r ' ll/y =(M. 00 A) 2 *F|V<H 

• Surface conditions 

• Blades (B): V<J>1T*0 

• Wakes (W):[V$]-n s O,ll<I>J = constant 

• Shocks (Sh):8>V<t«]| n* |4>1*0 

• Far-field conditions 

• Uniform flow conditions 

• Analytic solutions for steady disturbances 
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THE LINEARIZED UNSTEADY FLOW - 1 


The field equations governing the first-order unsteady perturbation of a nonlinear homentropic 
and irrotational steady flow are determined from the full time-dependent mass, momentum and 
entropy transport equations and the thermodynamic equation relating the entropy, pressure and 
density of a perfect gas. These equations can be cast in a convenient form by introducing the 
Goldstein (ref. 2) velocity decomposition, i.e., V = v R + S/<f>. The rotational component of the 
unsteady velocity (u H ) is divergence-free far upstream of the blade row and it is independent of 
the pressure fluctuation (p); the irrotational component (S/(f>) * s related directly to the unsteady 
pressure fluctuation. The resulting field equations for the first-order entropy (.s), rotational velocity 
(v R ) and velocity potential (cf>) are given below. First-order partial differential equations describe 
the transport of entropy and rotational velocity through the blade row. The unsteady potential 
is governed by a second-order equation, which is locally elliptic at field points at which the local 
steady Mach number (M) is less than one and locally hyperbolic at those points at which M > 1. 
Note that the rotational velocity provides a forcing function term to the potential equation. Also, 
if there are no entropy and rotational velocity fluctuations at inlet, then only a single field equation 
must be solved to determine the unsteady flow. 

• Velocity decomposition .• 7 = 7 R ♦ v<#> 

where T R = f (s-^v^VS)) and p=-p ^ 


• Field equations 



^-(v r -sv<!>/ 2) ♦ [(v R - sv<J>/2)-v] v4> a 0 
Dt (A " 2 of >-?"'*• •£v R) 
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THE LINEARIZED UNSTEADY FLOW - 2 


Conditions on the linearized unsteady perturbation at Blade ( B ), wake (IT) and shock (Sh) 
mean positions are obtained in a similar fashion, i.e., by substituting the asymptotic, Taylor and 
surface-vector series expansions into the full time-dependent surface conditions and equating terms 
of like order in e. The resulting first-order flow tangency, wake-jump (continuity of normal velocity 
and pressure) and shock-jump (conservation of mass and tangential momentum) conditions are 
indicated schematically below. Note that the blade displacement (rg) is prescribed but that the 
normal(to the shock) component of the shock displacement (r*sh • n) must be determined as part 
of the unsteady solution. Wake displacements have no 'impact on the linearized unsteady problem. 
In addition to the surface conditions, we require information on the unsteady flows far upstream 
and far downstream from the blade row. In these regions the linearized unsteady equations reduce 
to the constant coefficient equations of classical linearized theory, and analytic solutions for the 
velocity potential fluctuations (cf> p ) due to acoustic response disturbances and the far-downstream 
potential fluctuations (<f> R ) associated with the vortical or rotational velocity disturbances can be 
determined. These analytic far-field solutions can be matched to near-field numerical solutions, and 
they thereby serve to complete the specification of the linearized unsteady boundary-value problem 
(ref. 3). 


• Surface conditions 

• Blades: <f>„ * f (4>,7 B ,v R ) 

•Wakes: ] * [pJ =0 


• Shocks: |[F(<J>)«#> + G (*)<#>„] = f (*,7^ n .v R ) 

•sh,n * ” 1^1/ IM, I s ]l * [v R J =0 

• Far -field conditions 

s -co »^-a> Prescribed 


* 




*,N 


+ +0UT + * 


Prescribed 
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NUMERICAL SOLUTION DOMAIN 


The foregoing steady and linearized unsteady boundary-value problems account for the effects 
of blade geometry, mean blade loading and transonic, including moving shock, phenomena on the 
unsteady fluctuations arising from small-amplitude harmonic excitations. The unsteady equations 
are linear, time-independent and contain variable coefficients which depend on a fully nonlinear ho- 
mentropic and irrotational steady background flow. Numerical resolutions of the nonlinear steady 
and the linearized unsteady problems are required to determine the aerodynamic response informa- 
tion needed for aeroelastic and aeroacoustic applications. Because of the cascade geometry and the 
assumed form of the unsteady excitations (i.e., periodic in t and 77), such resolutions are required 
only over a single extended blade-passage region. In addition, since analytic far-field unsteady so- 
lutions have been determined, the numerical solution domain can be restricted further to a single 
extended blade-passage region of finite extent in the axial direction as shown below. 




W: Continuity of pressure and normal 


velocity across wakes 


Sh: Conservation of mass and tangential 


momentum across shocks 


P: Periodicity upstream 
M: Near/far field matching 
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AERODYNAMIC RESPONSE AT A BLADE SURFACE 


For aeroelastic and aeroacoustic applications, solutions to the nonlinear steady and the linearized 
unsteady boundary-value problems are required to predict the unsteady aerodynamic response at 
a moving blade surface (i.e., the unsteady surface pressures and the global unsteady airloads) and 
in the far-field (i.e., the unsteady pressure fluctuations), respectively. In particular, the pressure 
acting at the instantaneous position (5) of a given blade surface is made up of two components: a 
harmonic component (px) which is determined by the steady ($) and the linearized unsteady (<f>) 
potentials and the prescribed blade displacement (fa), and an anharmonic component (p Sh ) which 
is caused by the motion (rsh ,b) of a shock along the blade surface. The anharmonic pressure is 
determined by analytically continuing the solution to the steady boundary- value problem from the 
mean to the instantaneous shock location (ref. 4). Although the pressure disturbance ps is not 
everywhere harmonic, its regions of an harmonicity are small. Consequently, the first order global 
unsteady airloads are harmonic in time (ref. 5). In particular, if each two-dimensional blade section 
undergoes a pitching oscillation (re = a x R v ) about an axis fixed to the blade, the first-harmonic 
unsteady moment is determined by integrating the product of the first-harmonic component of the 
unsteady surface pressure and R p • f over the mean blade surface and subtracting a term consisting 
of the product of the steady pressure jump across the shock, R v • r and the shock displacement 
along the blade surface. 


Surface pressure 

P^= P H ($, 0. n> + P Sh <$- r Sh, B’tH - 



Unsteady moment 

m = /PHRp* d ^- r sh,B I t P B^P ,ir )lB + - 
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NUMERICAL RESULTS 


At this point we have completed our description of the unsteady aerodynamic model and pro- 
ceed to present numerical results to partially illustrate the status of numerical procedures for solving 
the nonlinear steady and linearized unsteady boundary-value problems and to demonstrate several 
important effects associated with nonuniform steady flow on the aerodynamic response at a moving 
blade surface. We refer the reader to refs. (6 and 7) for a description of the numerical procedures 
used. We will present results for two-dimensional compressor- or fan-type cascades operating at 
subsonic inlet and exit conditions. Theoretical results for steady surface Mach number (M) distri- 
butions, first-harmonic unsteady pressure-difference (A ph) distributions and unsteady aerodynamic 
moments (m) will be presented for blades undergoing pure pitching (torsional) motions with a = 1 , 0 
about their midchords. The stability of such motions depends upon the sign of the out-of-phase 
moment (m/). If m/ > 0, the airstream supplies energy to the blade motion, and this motion is 
unstable according to linearized theory. We will consider a subsonic cascade of NACA 0012 air- 
foils to illustrate the effects of a relatively thick, blunt-nosed blade geometry and variable mean 
incidence on the unsteady response and a subsonic/transonic cascade of 5% thick flat-bottomed 
double-circular-arc (DCA) airfoils to illustrate the effects of mean blade loading and transonic phe- 
nomena on the response at high subsonic inlet Mach number. For purposes of comparison results 
for flat-plate cascades, operating in uniform mean flows will be included along with those for the 
NACA 0012 and DCA cascades. The example cascades each have a stagger angle 0 of 45 deg and 
a unit gap/chord ratio ( G = 1). The steady flows through the NACA 0012 and DCA cascades have 
been determined by imposing a zero-load condition at sharp blade edges. 


• Torsional vibrations about midchord: r B =(ce x Rp) 

• Response parameters: Ap H and m 

• Subsonic flow: NACA 0012 airfoils 

• Subsonic/transonic flow: DCA airfoils 

mu** 


619 



SUBSONIC FLOW: EFFECT OF INCIDENCE - 1 


Steady and unsteady flows through the staggered cascade of modified NACA 0012 airfoils have 
been determined for an inlet Mach number (M_ oc ) of 0.6 and four inlet flow angles. The predicted 
surface Mach number distributions for fl-oo = 48, 50, 52 and 54 deg are shown below. The calculated 
exit Mach numbers are respectively 0.595, 0.557, 0.522 and 0.490, and in each case the calculated 
exit flow angle is approximately 47.7 deg. These steady flows are entirely subsonic with a peak Mach 
number of 0.789 occurring at x = 0.113 on the pressure (lower) surface of the blade for = 48 
deg, and 0.8, 0.86 and 0.96 occurring at x = 0.07, 0.05 and 0.03 on the suction (upper) surface for 
Cl-oo = 50, 52 and 54 deg, respectively. In each case the mean flow stagnates within 0.2% of blade 
chord downstream from the leading edge. 
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SUBSONIC FLOW: EFFECT OF INCIDENCE - 2 


Unsteady response predictions for the NACA 0012 cascade and for a flat-plate cascade operating 
at = 0.6 and fl _ <*, = 45° are shown on this and the next twofigures. Shown beloware unsteady 

pressure-difference distributions and aerodynamic moments for the reference (m = 0) NACA 0012 
and flat-plate blades undergoing unit-frequency pitching motions at a = 90 deg. The unsteady 
pressure difference is singular and behaves like a multiple of x ~ near the leading edge of the flat- 
plate airfoil. In contrast, the unsteady pressure is analytic in the vicinity of the rounded leading 
edge of the NACA 0012 blade. In this case both the real and imaginary components of the unsteady 
pressure difference are zero at the leading edge and reach local extrema very close to the leading 
edge. The results indicate that the coupling between the steady and unsteady flows, due to blade 
geometry and mean loading, leads to a reduction in the out-of-phase pressure difference, Im{Ap(x)}, 
over a forward part of the NACA 0012 blade and, therefore, a reduction in the out-of-phase moment 
opposing the blade motion. 


Unsteady pressure difference distributions: NACA 0012 cascade; 

0 = 450, G = 1.0, M _ 00 = 0.6, co = 1.0, o = 90 °. 

NACA 0012 cascade, flat-plate cascade 
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SUBSONIC FLOW: EFFECT OF INCIDENCE - 3 


Unsteady moments acting on the reference blade of the NACA 0012 cascade operating at inlet 
flow angles of 48, 50, 52 and 54 deg and on the reference blade of the flat-plate cascade operating 
at Q-oo = 45 deg are shown below for blades undergoing unit-frequency torsional vibrations over 
the entire range of interblade phase angles, i.e., cre[— 7r,7r]. The abrupt changes in the moment 
curves are indicative of an acoustic resonance. The blade motions are superresonant (i.e., acoustic 
response disturbances persist in the far field and carry energy away from the blade row) at interblade 
phase angles lying between the lowest and highest resonant phase angles and subresonant (acoustic 
response disturbances attenuate with increasing distance from the blade row) at the interblade phase 
angles below the lowest and above the highest resonant phase angles. The blade motions considered 
below are stable but the NACA 0012 results indicate that the effect of mean blade loading tends 
to be destabilizing. Note that for a given a the out-of-phase moment moves closer to the stability 
boundary as the inlet flow angle is increased. 


Unsteady moment vs. interblade phase angle: 0 = 45°, G = 1.0, M_oo = 0.6. 
NACA 0012 cascade, flat-plate cascade 
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SUBSONIC FLOW: EFFECT OF FREQUENCY 


The effect of frequency on the out-of-phase component of the unsteady moment due to torsion 
about midchord is illustrated below for the NACA 0012 cascade operating at oo = 54 deg and 
for the flat-plate cascade operating at ft_oo = 45 deg. The NACA 0012 blades experience a region 
of subresonant torsional instability for u? = 0.25, 0.5 and 0.75, with the extent of this region 
decreasing with increasing frequency. The subresonant torsional motions of the flat-plate cascade 
are unstable only at the lowest frequency considered, i.e., u = 0.25. Thus the nonuniform flow 
through the NACA 0012 cascade extends the frequency range over which the blades are susceptible 
to a torsional instability. 


Unsteady moment vs. interblade phase angle: 
0 = 45°, G = 1.0, M _ oo = 0.6. 


NACA 0012 cascade, Q_oo = 54° Flat-plate cascade, Q_oo = 45° 



623 



SUBSONIC/TRANSONIC FLOW: EFFECT OF MACH NUMBER - 1 


We now consider the staggered cascade of sharp-edged double-circular-arc airfoils. In particu- 
lar, the airfoils have flat lower surfaces, circular arc upper surfaces and maximum thicknesses at 
midchord of 0.05. Full-potential steady and linearized unsteady flows through this example config- 
uration have been determined for inlet Mach numbers of 0.5, 0.7, 0.8 and 0.9. The steady flows 
have been determined by imposing a zero-load requirement at blade leading and trailing edges. As 
a consequence, only the inlet Mach number is prescribed with the remaining inlet and exit param- 
eters determined as part of the steady flow solution. Numerical results for this configuration along 
with those for a corresponding flat-plate cascade (0 = 45 deg, G = 1) are given in this and in 
the following three figures. Shown below are the predicted surface Mach number distributions for 
the example DCA cascade. For the prescribed inlet Mach numbers stated above, the calculated 
exit Mach numbers are 0.43, 0.57, 0.62 and 0.64, respectively. In addition, the calculated inlet flow 
angles are 49.0, 49.2, 49.4, and 49.6 deg, respectively, and in each case the calculated exit flow angle 
is approximately 43.0 deg. The steady flows at A/_oo = 0.5, A/_ <» = 0.7 and A/_oo = 0.8 are entirely 
subsonic with the maximum suction-surface Mach numbers of 0.561, 0.804 and 0.941 occurring at, 
respectively, 40.8, 38.5 and 36.5% of blade chord downstream from the leading edge. The steady 
flow at A/_ oc = 0.9 is transonic with the supersonic region extending from 18.5 to 52.5% of blade 
chord along the suction surface and terminating at a shock discontinuity. The Mach numbers at 
the foot of the shock are 1.193 on the upstream or supersonic side and 0.871 on the downstream or 
subsonic side. 


Surface Mach number distributions: DCA cascade; 0 = 45°, G = 1.0 
Suction surface Pressure surface 
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SUBSONIC/TRANSONIC FLOW: EFFECT OF MACH NUMBER - 2 

The effect of Mach number on the response to in-phase (<r = 0 deg) unit-frequency torsional blade 
vibrations is illustrated below for the DCA and flat-plate cascades. The blade motions at a — Odeg 
are superresonant. For the flows at M - = 0.5, 0.7 and 0.8 two acoustic waves persist in the far 
field — one upstream and one downstream— and propagate away from the blade row. For the DCA 
cascade operating at M_ <» = 0.9 there are three such waves — two upstream and one downstream. 
Finally, for the flat-plate cascade operating at M_ = 0.9 there are four such waves — two upstream 
and two downstream. The out-of-phase pressure-difference distributions and unsteady moments for 
the reference DCA and flat-plate blades reflect this change in character of the acoustic response in 
the far field, since the trends indicated by the results for M_oo = 0.5, 0.7 and 0.8 are not sustained 
at M-oo = 0.9. Also, a comparison of the DCA and flat-plate pressure-difference curves for in-phase 
motions suggests that the influence of mean flow gradients on the unsteady aerodynamic response 
becomes more pronounced with increasing Mach number. The pressure difference distributions 
for the DCA and flat-plate blades are very similar for the two lower inlet Mach number, differ 
somewhat for M _<*, = 0.8 and differ substantially for M_ «*, = 0.9. The differences at M_ = 0.8 
can be attributed to the relatively large gradients in the subsonic mean flow that occur along the 
suction surface of each DCA blade. The substantial differences at M_oo = 0.9 are caused by the 
transonic effects associated with the DCA cascade and by the different far-field acoustic response 
environments produced by the two cascades. 


Unsteady pressure-difference distributions and unsteady moments: 

0 = 45°, G = 1.0, CO = 1.0, 0 = 0°. 
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SUBSONIC/TRANSONIC FLOW: EFFECT OF FREQUENCY 


Shown below are first-harmonic, out-of-phase, pressure-difference distributions and aerodynamic 
moments for the example DCA and flat-plate cascades operating at an inlet Mach number of 0.9. 
Here the blades are undergoing out-of-phase (<r = 180 deg) torsional vibrations about midchord 
at different prescribed frequencies. Recall that for a discontinuous transonic flow there are two 
contributions to the first-harmonic unsteady moment: one arising from the harmonic unsteady 
surface-pressure response and the other from the anharmonic surface pressures produced by shock 
motion. However, for the DCA cascade at M_ 0 0 = 0.9, the mean shock location is only slightly aft 
of blade midchord and, therefore, the anharmonic surface pressures make only a small contribution 
to the unsteady moment. A comparison of the DCA and flat-plate results depicted below indicates 
the dramatic impact of transonic mean-flow phenomena on unsteady aerodynamic response. A 
second interesting feature indicated by these results is the change in the unsteady moment behavior 
as the blade vibration frequency is increased from 0.25 to 0.5. This change occurs because the 
out-of-phase blade motions of the DCA and flat-plate cascades are subresonant for u> = 0.1 and 
0.25 and superresonant for u> = 0.5, 0.75 and 1.0, and this change in the far-field acoustic response 
has an important impact on the unsteady aerodynamic response at a blade surface. 


Unsteady pressure-difference distributions: 
0 = 45°, G = 1.0, M_oo = 0.9, O = 180°. 



0 0.2 0.4 0.6 0.8 1.0 
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SUBSONIC/TRANSONIC FLOW 


Finally, predictions of the aerodynamic moment versus interblade phase angle are shown for 
unit frequency torsional blade vibrations for the flat-plate and DCA cascades. Here M_ ^ = 0.9 
and those angles at which an acoustic resonance occurs are indicated by the arrows at the top of 
each figure. The unit-frequency torsional blade motions of the two cascades are stable (i.e., m/ < 0 
for all a ), but the behaviors of the flat-plate and DCA moment responses vs. interblade phase 
angle are quite different. These differences occur not only because the mean flow through the flat- 
plate cascade is entirely subsonic while that through the DCA cascade is transonic with a shock 
discontinuity, but also because of the substantial difference between the exit Mach numbers for the 
flat-plate (M 0 0 = 0.9) and the DCA (M 0 0 = 0.64) cascades. This difference implies that the two 
cascades operate in very different far-downstream acoustic response environments over almost the 
entire range of interblade phase angles. 


Unsteady moment vs. interblade phase angle; m R , m, 

0 = 45°, G = 1 , M_ oo = 0.9, cu = 1.0. 



o (deg) o (deg) 
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SUMMARY 


The linearized unsteady aerodynamic theory outlined above accounts for the effects of real blade 
geometry, mean blade loading and operation at transonic Mach numbers on the unsteady aerody- 
namic response produced by the blades of an isolated two-dimensional cascade. This theory has 
been developed to meet the requirements of turbomachinery aeroelastic designers, but it should 
also be useful for aeroacoustic design applications. The unsteady flow is regarded as a small pertur- 
bation of a fully nonuniform isentropic and irrotational mean or steady flow, which is produced by 
small-amplitude temporally and spatially (in the cascade direction) periodic structural (blade mo- 
tions) and external aerodynamic (incident entropic, vortical and acoustic disturbances) excitations. 
Thus the steady flow is determined as a solution of a full-potential boundary-value problem and the 
linearized unsteady flow as a solution of a time-independent, linear, variable-coefficient, boundary- 
value problem in which the variable coefficients depend on the underlying mean or steady flow. 

Response predictions have been presented for the blades of compressor- and fan-type cascades 
undergoing pure torsional motions. In these examples there are no incident entropy or rotational 
velocity fluctuations and, therefore, only a single field equation must be solved to determine the 
linearized unsteady flow field. The numerical results demonstrate, to some extent, the status of nu- 
merical field methods for solving the nonlinear steady and the linear, variable-coefficient, unsteady, 
boundary-value problems and illustrate partially the effects of blade geometry, mean incidence, 
shock phenomena and differences between inlet and exit free-stream conditions on the unsteady 
response at blade surfaces. 


• Linearized unsteady aerodynamic analysis 
e Effects of: 

• Blade geometry 

• Blade loading 

• Shocks and their motions 

• High frequency unsteady motions 

e Blade flutter prediction 
e Subsonic /transonic Mach numbers 
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FUTURE DIRECTIONS 


Linearizations relative to nonuniform steady flows offer great potential for meeting the needs 
of aeroelastic (or aeroacoustic) designers for efficient unsteady aerodynamic analyses that contain 
much of the essential physics associated with turbomachinery flow fields. However, before this po- 
tential can be fully realized, significant improvements in numerical solution methods for both the 
steady and linearized unsteady flows must be achieved so that reliable response information can 
be provided over the wide range of geometric configurations and flow conditions at which blade 
vibrations are of practical concern. In particular, unsteady aerodynamic analyses intended for tur- 
bomachinery aeroelastic predictions must be applicable to fan, compressor apd turbine cascades, 
low subsonic through low supersonic Mach number operation and moderate through high frequency 
structural and external aerodynamic excitations. Some needed capabilities include the ability to 
predict transonic flows (i.e., subsonic flows with imbedded supersonic regions) through fan and 
compressor cascades operating at high positive or negative mean incidence, supersonic flows with 
complicated moving shock patterns and the high frequency unsteady flows caused by incident exter- 
nal aerodynamic disturbances. Major advances in our ability to predict turbomachinery aeroelastic 
and aeroacoustic behavior should result if future research is directed toward including the effects 
of strong viscid/inviscid interactions and possibly largje-scale flow-separations within a linearized 
unsteady aerodynamic framework. Ultimately, linearized analyses which account for nonuniform 
steady flow and viscid/inviscid interaction phenomena must be extended to treat three-dimensional 
flows. 


• Subsonic/transonic flows at high incidence 


• Forced aerodynamic excitations 

•^>=f(v?s) 

• Supersonic Mach numbers 

• Viscous separation phenomena 
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TRANSONIC AERODYNAMICS AND FLUTTER 


In the past decade there has been much activity in the development of computational methods for 
the analysis of unsteady transonic aerodynamics about airfoils and wings. The upper left figure 
illustrates significant features which must be addressed in the treatment of computational 
transonic unsteady aerodynamics. On the plot of equivalent airspeed versus Mach number, lines 
of constant altitude are straight lines through the origin with decreasing altitudes represented 
by steeper slopes. The flight envelope, typically set by the maximum limit speed and a typical 
flutter boundary curve, characterized by the flutter speed gradually dropping to a minimum in 
the transonic speed range followed by a rapid upward rise, is shown. The ability to predict this 
minimum, termed the transonic flutter dip, is of great importance in design, since the flutter 
boundary must be shown by a combination of analysis and flight test to be outside the flight 
envelope by a margin of at least 15 percent in equivalent airspeed for military aircraft. 

The upper right figure indicates the flow regions for an aircraft on a plot of lift coefficient 
versus Mach number. Flows which are predominantly attached or separated are designated as 
type I and III respectively, while mixed attached and separated flows are designated type II. For 
aeroelastic problems the boundary of the type II flows will be enlarged over that for steady 
flows since a vibrating airfoil or wing may exhibit alternating attached and separated flow for 
sensitive conditions. The "picket fence" in the mixed flow region has been added to emphasize 
the possibility of "nonclassical" aeroelastic effects in this region. 

The diagram in the lower left of the figure illustrates the sequence of events occuring in air 
combat maneuvers. Upon the decision to engage, a maneuver is initiated with the objective of 
achieving maximum turn rate. This leads, in turn, to pull-up and turn at the structural limit 
load, decelerating at limit load to the intersection with the maximum lift coefficient curve,' 
holding this "corner" condition until the pointing objective is achieved and completion of 
engagement and pull-out occurs. These maneuvers, encompassing the complete fight envelope, 
involve rapid transitions between type I, II, and III flow conditions. 

Further features of transonic flutter are illustrated in the lower right diagram. Dynamic 
pressure at flutter tends to decrease with increasing Mach number to a minimum "critical 
flutter point" value in the transonic speed range. At subsonic speeds the flow can be reasonably 
assumed to be attached (type I) at flutter and linear theory is well calibrated for flutter 
analysis. At transonic speeds the situation is complicated by the onset of flow separation (type 
II flow) and linear theory must be used with caution. The low damping region indicated in the 
figure indicates the potential for nonclassical aeroelastic response and instabilities which may 
be encountered. 

References : 

Edwards, J. W.; and Thomas, J. L.: Computational Methods for Unsteady Transonic Flows, AIAA 
Paper No. 87-0107. 
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TRANSONIC AERODYNAMICS AND FLUTTER 


GRAPHICAL REPRESENTATION OF MINIMUM 
REQUIRED FLUTTER MARGIN 


CHARACTERISTICS OF ATTACHED AND 
SEPARATED FLOW FOR COMPLETE AIRCRAFT 
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COMPUTATIONAL AEROELASTICITY CHALLENGES 


This figure illustrates several types of aeroelastic response which have been encountered and 
which offer challenges for computational methods. The four cases illustrate problem areas 
encountered near the boundaries of aircraft flight envelopes, as operating conditions change 
from high speed, low angle conditions to lower speed, higher angle conditions. The nonclassical 
aeroelastic response observed on the DAST ARW-2 wing model (upper left) is a region of high 
dynamic response at nearly constant Mach number which was encountered at dynamic pressures 
well below those for which flutter was predicted. The motion is of the limit-amplitude type and 
the response is believed to be associated with flow separation and reattachment over the 
supercritical wing (type II flow). 

The upper right figure illustrates wing/store limited amplitude oscillations experienced by 
modern, high performance aircraft under various loading and maneuvering conditions at 
transonic Mach numbers. Such oscillations can result in limitations on vehicle performance. 
The conditions for which this response occurs appear to be near the onset of type II mixed flow. 
The response typically increases for maneuvering flight conditions. 

Dynamic vortex-structure interactions causing wing oscillations have been observed on a 
bomber type aircraft for high wing sweep conditions during wind-up turn maneuvers (lower 
left). The flow involves the interaction of the wing vortex system with the first wing bending 
mode and occurs over a wide Mach number range (0.6 - 0.95) at angles of attack of 7 - 9 
degrees. 

At higher angles, interaction of forebody and wing vortex systems with aft vehicle components 
results in vortex-induced buffet loads, illustrated in the lower right figure. The figure shows 
the operating conditions for which tail buffet may occur on a high performance fighter. Buffet 
of horizontal tails can occur at intermediate angles of attack and is a result of the vortex system 
encountering the horizontal tail lifting surface. As angle of attack increases, the location of 
vortex bursting moves upstream in the wake. Loss of lift is associated with the burst location 
reaching the vicinity of the aircraft, and vertical tail surfaces located in such regions can 
experience severe dynamic loads. 
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COMPUTATIONAL AEROELASTICITY CHALLENGES 


NOVEL SHOCK-INDUCED INSTABILITIES WING/STORE LIMITED AMPLITUDE FLUTTER 



Linear theory ftuttar boundary 
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COMPUTER RESOURCE REQUIREMENTS FOR FLUTTER ANALYSIS 


This table indicates the computer resources required to perform a flutter analysis of a complete 
aircraft configuration at one Mach number. Time-marching transient aeroelastic response 
calculations are used to determine the flutter condition. This involves, on average, four 
response calculations: two to calculate steady flow field conditions and two transient responses 
bracketing the flutter speed. Modal frequency and damping estimates from the responses are 
determined and the flutter speed interpolated from the damping estimates. Calculations have 
been performed for a complete aircraft configuration with a transonic small disturbance (TSD) 
potential code using 750,000 grid points. The calculation of one flutter point for this case on 
the CDC VPS-32 computer would require 2.3 CPU hours. Estimates of similar calculations 
using the full Navier-Stokes equations would require 77.8 CPU hours. Conditions for this 
estimate are a Reynolds number of 10 million, 7 million grid points and an assumed 
computational speed of 100 million floating point operations per second (MFLOPS). 

References : 

Whitlow, Woodrow, Jr.: Computational Unsteady Aerodynamics for Aeroelastic Analysis, NASA 
TM 100523, December 1987. 


COMPUTER RESOURCE REQUIREMENTS TO DETERMINE FLUTTER POINT 
AT A SPECIFIED MACH NUMBER 

(4000 TIME STEPS PER FLUTTER POINT) 





CPU HOURS 
(YPS-32) 

COMPLETE AIRCRAFT 

TSD 

0.75M 

2.3* 

COMPLETE AIRCRAFT 

FULL NAVIER-STOKES 

7.00M 

77.8“ 


(RE = 10 MILLION) 




‘BASED ON ACTUAL CASES 
“ASSUMES COMPUTATIONAL SPEED OF 100 MFLOPS 
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COMPUTER RESOURCE REQUIREMENTS FOR 
COMPLETE FLUTTER BOUNDARY 


This table summarizes computational requirements for flutter calculations of a 
wing/body/canard configuration on the CDC VPS-32 computer operating at 100 MFLOPS and on 
the NAS CRAY II computer operating at 250 MFLOPS. Again, four response calculations per 
flutter point are assumed. It is assumed that ten flutter points will be calculated to define the 
flutter boundary versus Mach number. The left hand column indicates the difficulty of the 
flowfield calculation as defined in figure 1; type I for attached flows, type II for mixed 
(alternately separated and attached) flows and type III for fully separated flows. The second 
column indicates the fluid dynamic equation level needed to accurately model the flow physics of 
the problem. Note that two-dimensional strip boundary layer models are assumed for 
interactive viscous-inviscid calculation's for the potential and Euler equation methods. It is 
anticipated that potential equation models will be adequate for flutter calculations of type I 
attached flow conditions and may also be quite useful for some type II mixed flow cases. Full 
potential equation codes will require about 50 percent more computer resources than TSD 
methods due to the necessity of conforming, moving grids, among other considerations. Euler 
equation methods should also be adequate for these conditions and, in addition, be able to treat 
more difficult type III fully separated flows. Euler equation methods are estimated to require 
approximately twice the resources of TSD methods. The full Navier-Stokes equations, which 
should only be required for type II and III flows require approximately 30 times the resources 
of the Euler equations (at a Reynolds number of 1 00 million). 


WING/BODY/CANARD CONFIGURATION 
10 MACH NUMBERS (40 CASES) PER ANALYSIS 

OPS OPS 

TIME = (GRID FTS)X— ) 


FLOW REGION 

FLOW MODEL 

VPS-32 (100 MFLOPS) 

NAS (250 MFLOPS) 

1, MAYBE II 

TSD WITH 2-D 
STRIP BOUNDARY LAYER 

30 HOURS 

12 HOURS 

1, MAYBE II 

POTENTIAL WITH 2-D 
STRIP BOUNDARY LAYER 

45 HOURS 

18 HOURS 

1, II, MAYBE III 

EULER WITH 2-D 
STRIP BOUNDARY LAYER 

65 HOURS 

26 HOURS 

II, III 

NAVIER-STOKES 
(RE = 108) 

1611 HOURS 

644 HOURS 
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